From d04668198632a57639b972cfaddbffe167243c21 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Wed, 26 Aug 2026 18:08:37 +0200 Subject: [PATCH 1/9] added bbse hard and soft and classifier adaptation wrapper --- TODO.txt | 4 +- docs/source/manuals.rst | 1 + docs/source/manuals/methods.md | 54 ++++ docs/source/quapy.classification.rst | 8 + quapy/classification/__init__.py | 1 + quapy/functional.py | 105 ++----- quapy/method/__init__.py | 6 + quapy/method/_liep_draft.py | 261 ++++++++++++++++ quapy/method/aggregative.py | 433 +++++++++++++++++++++++++-- quapy/tests/test_methods.py | 26 +- 10 files changed, 791 insertions(+), 108 deletions(-) create mode 100644 quapy/method/_liep_draft.py diff --git a/TODO.txt b/TODO.txt index 17604ab..12f1f59 100644 --- a/TODO.txt +++ b/TODO.txt @@ -19,4 +19,6 @@ scale each value by per-class thresholds, i.e., [0.33*0.1, 0.33*1, 0.33*1]/sum. - [TODO] add Friedman's method and DeBias - [TODO] check ignore warning stuff check https://docs.python.org/3/library/warnings.html#temporarily-suppressing-warnings -- [TODO] nmd and md are not selectable from qp.evaluation.evaluate as a string \ No newline at end of file +- [TODO] nmd and md are not selectable from qp.evaluation.evaluate as a string +- [TODO] add https://www.kaggle.com/datasets/Cornell-University/arxiv dataset? in the paper "Online Adaptation to + Label Distribution Shift" they say it provides a natural label shift over time \ No newline at end of file diff --git a/docs/source/manuals.rst b/docs/source/manuals.rst index 97e0628..f81ab7c 100644 --- a/docs/source/manuals.rst +++ b/docs/source/manuals.rst @@ -7,6 +7,7 @@ Manuals manuals/datasets manuals/evaluation + manuals/label-shift-adaptation manuals/methods manuals/model-selection manuals/plotting diff --git a/docs/source/manuals/methods.md b/docs/source/manuals/methods.md index b721957..5be168f 100644 --- a/docs/source/manuals/methods.md +++ b/docs/source/manuals/methods.md @@ -323,6 +323,60 @@ model.fit(*train.Xy) estim_prevalence = model.predict(test.X) ``` +RLLS computes its importance weights directly (it implements the `ImportanceWeightQuantifier` interface); see +the {ref}`Label Shift Adaptation manual ` for how to +access these weights, or use them to adapt a classifier itself rather than only estimating prevalence. + +### Black Box Shift Estimation (BBSE) + +`BBSEhard` and `BBSEsoft` are available at `qp.method.aggregative.BBSEhard` and +`qp.method.aggregative.BBSEsoft`, respectively, and implement the Black Box Shift Estimator +proposed in: + +_Lipton, Z., Wang, Y. X., & Smola, A. (2018, July). Detecting and correcting for label shift +with black box predictors. In International conference on machine learning +(pp. 3122-3130). PMLR._ ([link to paper](https://proceedings.mlr.press/v80/lipton18a.html)) + +BBSE is similar in spirit to ACC and PACC in that it exploits the label-shift invariance of +`P(hat{Y}|Y)` to correct for the change in class prevalence between the training and the test +distributions. However, while ACC solves the linear system `q = Mp` (with `M` the matrix of +class-conditional misclassification rates and `p` the sought prevalence vector), BBSE instead +solves `q = Cw` for the importance-weight vector `w`, with `w_i = Q(i)/P(i)` the ratio between +the target and the source class priors, and `C` the joint-distribution matrix +`C_ij = P(hat{Y}=i, Y=j)` estimated on a validation split. The target prevalence estimate is +then recovered as `Q(y) = w_y * P(y)`. + +`BBSEhard` estimates `C` from crisp classifier predictions (i.e., a standard confusion matrix, +normalized to sum to 1), while `BBSEsoft` estimates it from the classifier's posterior +probabilities instead, in the same spirit in which PACC generalizes ACC. + +```python +import quapy as qp +from quapy.method.aggregative import BBSEhard, BBSEsoft +from sklearn.linear_model import LogisticRegression + +train, test = qp.datasets.fetch_UCIBinaryDataset('haberman').train_test + +model = BBSEhard(LogisticRegression(max_iter=2000), val_split=5) +model.fit(*train.Xy) +estim_prevalence = model.predict(test.X) + +# or, using posterior probabilities instead of crisp counts: +model = BBSEsoft(LogisticRegression(max_iter=2000), val_split=5) +model.fit(*train.Xy) +estim_prevalence = model.predict(test.X) +``` + +As with ACC and RLLS, both variants require validation predictions and therefore expect +`val_split` to be set whenever `fit_classifier=True`. They also accept the same `solver` +(`"minimize"`, `"exact-raise"`, `"exact-cc"`) and `norm` (`"clip"`, `"mapsimplex"`, +`"condsoftmax"`) arguments discussed above for ACC/PACC. + +Like RLLS, both `BBSEhard` and `BBSEsoft` implement the `ImportanceWeightQuantifier` interface; see the +{ref}`Label Shift Adaptation manual ` for how to +access these weights directly, or use them to adapt a classifier itself rather than only estimating +prevalence. + ### Distribution Matching Distribution Matching (DM) methods search for the mixture parameter (the sought class prevalence values) diff --git a/docs/source/quapy.classification.rst b/docs/source/quapy.classification.rst index cfc7d9b..c07b72e 100644 --- a/docs/source/quapy.classification.rst +++ b/docs/source/quapy.classification.rst @@ -12,6 +12,14 @@ quapy.classification.calibration module :undoc-members: :show-inheritance: +quapy.classification.labelshift module +--------------------------------------- + +.. automodule:: quapy.classification.labelshift + :members: + :undoc-members: + :show-inheritance: + quapy.classification.methods module ----------------------------------- diff --git a/quapy/classification/__init__.py b/quapy/classification/__init__.py index b63e0eb..33551b7 100644 --- a/quapy/classification/__init__.py +++ b/quapy/classification/__init__.py @@ -1,3 +1,4 @@ from . import calibration +from . import labelshift from . import methods from . import svmperf diff --git a/quapy/functional.py b/quapy/functional.py index 392309e..84504ee 100644 --- a/quapy/functional.py +++ b/quapy/functional.py @@ -462,7 +462,7 @@ def argmin_prevalence(loss: Callable, raise NotImplementedError() -def optim_minimize(loss: Callable, n_classes: int, return_loss=False): +def optim_minimize(loss: Callable, n_classes: int, x0='uniform', bounds='simplex', constraints='simplex', return_loss=False): """ Searches for the optimal prevalence values, i.e., an `n_classes`-dimensional vector of the (`n_classes`-1)-simplex that yields the smallest lost. This optimization is carried out by means of a constrained search using scipy's @@ -470,19 +470,32 @@ def optim_minimize(loss: Callable, n_classes: int, return_loss=False): :param loss: (callable) the function to minimize :param n_classes: (int) the number of classes, i.e., the dimensionality of the prevalence vector + :param x0: initial solution; if the string 'uniform' is passed (default) then the initial solution is the + uniform distribution; otherwise, a valid object must be provided + :param bounds: the bounds of the search space; if the string `simplex` is passed (default) then the bounds + of a simplex of appropriate dimension is instantiated; otherwise, a valid tuple must be provided + :param constraints: the constraints of valid solutions; if the string `simplex` is passed (default) then the + constraint that a point must lie on the simplex is assumed; otherwise, valid constrains must be provided + (see scipy.optimize) :param return_loss: bool, if True, returns also the value of the loss (default is False). :return: (ndarray) the best prevalence vector found or a tuple which also contains the value of the loss if return_loss=True """ from scipy import optimize - # the initial point is set as the uniform distribution - uniform_distribution = uniform_prevalence(n_classes=n_classes) + if isinstance(x0, str) and x0=='uniform': + # the initial point is set as the uniform distribution + x0 = uniform_prevalence(n_classes=n_classes) - # solutions are bounded to those contained in the unit-simplex - bounds = tuple((0, 1) for _ in range(n_classes)) # values in [0,1] - constraints = ({'type': 'eq', 'fun': lambda x: 1 - sum(x)}) # values summing up to 1 - r = optimize.minimize(loss, x0=uniform_distribution, method='SLSQP', bounds=bounds, constraints=constraints) + if isinstance(bounds, str) and bounds=='simplex': + # solutions are bounded to those contained in the unit-simplex + bounds = tuple((0, 1) for _ in range(n_classes)) # values in [0,1] + + if isinstance(constraints, str) and constraints=='simplex': + # the point must lie on the simplex + constraints = ({'type': 'eq', 'fun': lambda x: 1 - sum(x)}) # values summing up to 1 + + r = optimize.minimize(loss, x0=x0, method='SLSQP', bounds=bounds, constraints=constraints) if return_loss: return r.x, r.fun @@ -623,84 +636,6 @@ def solve_adjustment_binary(prevalence_estim: ArrayLike, tpr: float, fpr: float, return adjusted -def solve_adjustment( - class_conditional_rates: np.ndarray, - unadjusted_counts: np.ndarray, - method: Literal["inversion", "invariant-ratio"], - solver: Literal["exact", "minimize", "exact-raise", "exact-cc"]) -> np.ndarray: - """ - Function that tries to solve for :math:`p` the equation :math:`q = M p`, where :math:`q` is the vector of - `unadjusted counts` (as estimated, e.g., via classify and count) with :math:`q_i` an estimate of - :math:`P(\hat{Y}=y_i)`, and where :math:`M` is the matrix of `class-conditional rates` with :math:`M_{ij}` an - estimate of :math:`P(\hat{Y}=y_i|Y=y_j)`. - - :param class_conditional_rates: array of shape `(n_classes, n_classes,)` with entry `(i,j)` being the estimate - of :math:`P(\hat{Y}=y_i|Y=y_j)`, that is, the probability that an instance that belongs to class :math:`y_j` - ends up being classified as belonging to class :math:`y_i` - - :param unadjusted_counts: array of shape `(n_classes,)` containing the unadjusted prevalence values (e.g., as - estimated by CC or PCC) - - :param str method: indicates the adjustment method to be used. Valid options are: - - * `inversion`: tries to solve the equation :math:`q = M p` as :math:`p = M^{-1} q` where - :math:`M^{-1}` is the matrix inversion of :math:`M`. This inversion may not exist in - degenerated cases. - * `invariant-ratio`: invariant ratio estimator of `Vaz et al. 2018 `_, - which replaces the last equation in :math:`M` with the normalization condition (i.e., that the sum of - all prevalence values must equal 1). - - :param str solver: the method to use for solving the system of linear equations. Valid options are: - - * `exact-raise`: tries to solve the system using matrix inversion. Raises an error if the matrix has rank - strictly lower than `n_classes`. - * `exact-cc`: if the matrix is not full rank, returns :math:`q` (i.e., the unadjusted counts) as the estimates - * `exact`: deprecated, defaults to 'exact-cc' (will be removed in future versions) - * `minimize`: minimizes a loss, so the solution always exists - """ - if solver == "exact": - warnings.warn( - "The 'exact' solver is deprecated. Use 'exact-raise' or 'exact-cc'", DeprecationWarning, stacklevel=2) - solver = "exact-cc" - - A = np.asarray(class_conditional_rates, dtype=float) - B = np.asarray(unadjusted_counts, dtype=float) - - if method == "inversion": - pass # We leave A and B unchanged - elif method == "invariant-ratio": - # Change the last equation to replace it with the normalization condition; - # copy first so this does not mutate the caller's arrays (np.asarray above - # returns the same object, not a copy, when the input is already float64) - A = A.copy() - B = B.copy() - A[-1, :] = 1.0 - B[-1] = 1.0 - else: - raise ValueError(f"unknown {method=}") - - if solver == "minimize": - def loss(prev): - return np.linalg.norm(A @ prev - B) - return optim_minimize(loss, n_classes=A.shape[0]) - elif solver in ["exact-raise", "exact-cc"]: - # Solvers based on matrix inversion, so we use try/except block - try: - return np.linalg.solve(A, B) - except np.linalg.LinAlgError: - # The matrix is not invertible. - # Depending on the solver, we either raise an error - # or return the classifier predictions without adjustment - if solver == "exact-raise": - raise - elif solver == "exact-cc": - return unadjusted_counts - else: - raise ValueError(f"Solver {solver} not known.") - else: - raise ValueError(f'unknown {solver=}') - - # ------------------------------------------------------------------------------------------ # Transformations from Compositional analysis # ------------------------------------------------------------------------------------------ diff --git a/quapy/method/__init__.py b/quapy/method/__init__.py index 0bcab12..695ac5f 100644 --- a/quapy/method/__init__.py +++ b/quapy/method/__init__.py @@ -15,6 +15,8 @@ AGGREGATIVE_METHODS = { aggregative.ACC, aggregative.PCC, aggregative.PACC, + aggregative.BBSEhard, + aggregative.BBSEsoft, aggregative.RLLS, aggregative.EMQ, aggregative.HDy, @@ -30,6 +32,7 @@ AGGREGATIVE_METHODS = { aggregative.KDEyML, aggregative.KDEyCS, aggregative.KDEyHD, + aggregative.LEIP, # aggregative.OneVsAllAggregative, confidence.BayesianCC, _bayesian.BayesianKDEy, @@ -54,12 +57,15 @@ MULTICLASS_METHODS = { aggregative.ACC, aggregative.PCC, aggregative.PACC, + aggregative.BBSEhard, + aggregative.BBSEsoft, aggregative.RLLS, aggregative.EMQ, aggregative.EDy, aggregative.KDEyML, aggregative.KDEyCS, aggregative.KDEyHD, + aggregative.LEIP, confidence.BayesianCC, _bayesian.BayesianKDEy, _bayesian.BayesianMAPLS, diff --git a/quapy/method/_liep_draft.py b/quapy/method/_liep_draft.py new file mode 100644 index 0000000..62c9bf7 --- /dev/null +++ b/quapy/method/_liep_draft.py @@ -0,0 +1,261 @@ +import numpy as np + + +def normalize(v, eps=1e-12): + v = np.asarray(v, dtype=float) + v = np.maximum(v, eps) + return v / v.sum() + + +def bayes_prior_update(probs, new_prior, source_prior, eps=1e-12): + """ + Applies the standard label-shift prior correction: + + p_new(y|x) ∝ p_old(y|x) * new_prior(y) / source_prior(y) + + Parameters + ---------- + probs : array, shape (n_samples, n_classes) or (n_classes,) + Probabilistic outputs of the classifier. + new_prior : array, shape (n_classes,) + Current estimate of the target prior. + source_prior : array, shape (n_classes,) + Source/training prior. + """ + probs = np.asarray(probs, dtype=float) + one_dim = probs.ndim == 1 + if one_dim: + probs = probs[None, :] + + new_prior = normalize(new_prior, eps=eps) + source_prior = normalize(source_prior, eps=eps) + + weights = new_prior / np.maximum(source_prior, eps) + updated = probs * weights[None, :] + updated = updated / np.maximum(updated.sum(axis=1, keepdims=True), eps) + + return updated[0] if one_dim else updated + + +def confusion_statistic(y_true, y_pred, n_classes, mode="recall", eps=1e-12): + """ + Computes the minimum diagonal statistic used to choose tau. + + mode="recall": + diag(C) / row sums, i.e. per-true-class recall. + + mode="precision": + diag(C) / column sums, i.e. per-predicted-class precision. + + The paper calls this quantity 'minimum recall', but also describes + a column-normalized confusion matrix, which would correspond more + closely to precision. I expose both options. + """ + y_true = np.asarray(y_true, dtype=int) + y_pred = np.asarray(y_pred, dtype=int) + + C = np.zeros((n_classes, n_classes), dtype=float) + for yt, yp in zip(y_true, y_pred): + C[yt, yp] += 1.0 + + if mode == "recall": + denom = C.sum(axis=1) + elif mode == "precision": + denom = C.sum(axis=0) + else: + raise ValueError("mode must be 'recall' or 'precision'.") + + diag = np.diag(C) + valid = denom > 0 + + if not np.any(valid): + raise ValueError("No valid classes found in the confusion matrix.") + + scores = diag[valid] / np.maximum(denom[valid], eps) + return float(np.min(scores)) + + +def choose_tau_from_validation( + target_probs, + val_probs, + y_val, + mode="recall", + retain_fraction=None, +): + """ + Chooses tau following the spirit of Section 4.1. + + The paper says tau is selected from the top n percentile of target + confidences, with n = min-recall * 100. To make the selected set A + contain approximately min-recall fraction of the target samples, we set: + + tau = quantile(max_probs, 1 - min_recall) + + so that about min_recall of the target points satisfy max_prob >= tau. + + If retain_fraction is provided, it overrides the validation-derived value. + """ + target_probs = np.asarray(target_probs, dtype=float) + val_probs = np.asarray(val_probs, dtype=float) + y_val = np.asarray(y_val, dtype=int) + + n_classes = target_probs.shape[1] + y_val_pred = val_probs.argmax(axis=1) + + if retain_fraction is None: + retain_fraction = confusion_statistic( + y_true=y_val, + y_pred=y_val_pred, + n_classes=n_classes, + mode=mode, + ) + + retain_fraction = float(np.clip(retain_fraction, 0.0, 1.0)) + + target_conf = target_probs.max(axis=1) + + if retain_fraction <= 0: + tau = np.inf + elif retain_fraction >= 1: + tau = -np.inf + else: + tau = np.quantile(target_conf, 1.0 - retain_fraction) + + return tau, retain_fraction + + +def leip( + target_probs, + source_prior, + tau=None, + val_probs=None, + y_val=None, + threshold_mode="recall", + count_smoothing=0.0, + eps=1e-12, + return_details=False, +): + """ + LEIP: Label shift Estimation with Incremental Prior update. + + Parameters + ---------- + target_probs : array, shape (n_target, n_classes) + Probabilistic classifier outputs on the target/test set. + source_prior : array, shape (n_classes,) + Source class prior p_s(y). + tau : float or None + Confidence threshold. If None, it is estimated using validation data. + val_probs : array, shape (n_val, n_classes), optional + Validation probabilistic outputs, required if tau is None. + y_val : array, shape (n_val,), optional + Validation labels, required if tau is None. + threshold_mode : {"recall", "precision"} + Statistic used to choose tau from validation data. + count_smoothing : float + Optional additive smoothing for pseudo-label counts. Set to 0.0 + for a closer implementation of the paper; use a small value such + as 1e-8 for extra numerical robustness. + eps : float + Numerical stabilizer. + return_details : bool + If True, returns diagnostic information. + + Returns + ------- + estimated_prior : array, shape (n_classes,) + Estimated target class distribution. + details : dict, optional + Returned only if return_details=True. + """ + target_probs = np.asarray(target_probs, dtype=float) + if target_probs.ndim != 2: + raise ValueError("target_probs must have shape (n_samples, n_classes).") + + n_target, n_classes = target_probs.shape + source_prior = normalize(source_prior, eps=eps) + + if n_classes != len(source_prior): + raise ValueError("source_prior must have one entry per class.") + + # Step 1: choose tau if needed + if tau is None: + if val_probs is None or y_val is None: + raise ValueError("val_probs and y_val are required when tau is None.") + + tau, retain_fraction = choose_tau_from_validation( + target_probs=target_probs, + val_probs=val_probs, + y_val=y_val, + mode=threshold_mode, + ) + else: + retain_fraction = None + + target_conf = target_probs.max(axis=1) + target_top = target_probs.argmax(axis=1) + + # Step 2: high-confidence set A + A_mask = target_conf >= tau + A_labels = target_top[A_mask] + + counts = np.full(n_classes, count_smoothing, dtype=float) + + if len(A_labels) > 0: + counts += np.bincount(A_labels, minlength=n_classes) + current_prior = counts / counts.sum() + else: + # Fallback if tau is too strict. + # One could also use classify-and-count over the full target set. + current_prior = source_prior.copy() + + # Step 3: low-confidence set B, sorted by decreasing confidence + B_indices = np.where(~A_mask)[0] + B_indices = B_indices[np.argsort(-target_conf[B_indices])] + + # Step 4: incremental pass over B + incremental_labels = [] + + for idx in B_indices: + corrected = bayes_prior_update( + probs=target_probs[idx], + new_prior=current_prior, + source_prior=source_prior, + eps=eps, + ) + pseudo_label = int(np.argmax(corrected)) + incremental_labels.append(pseudo_label) + + counts[pseudo_label] += 1.0 + current_prior = counts / counts.sum() + + estimated_intermediate_prior = current_prior.copy() + + # Step 5: final complete pass over all target instances + corrected_all = bayes_prior_update( + probs=target_probs, + new_prior=estimated_intermediate_prior, + source_prior=source_prior, + eps=eps, + ) + + final_labels = corrected_all.argmax(axis=1) + estimated_prior = np.bincount(final_labels, minlength=n_classes).astype(float) + estimated_prior /= estimated_prior.sum() + + if not return_details: + return estimated_prior + + details = { + "tau": tau, + "retain_fraction": retain_fraction, + "n_A": int(A_mask.sum()), + "n_B": int((~A_mask).sum()), + "A_mask": A_mask, + "intermediate_prior": estimated_intermediate_prior, + "final_labels": final_labels, + "corrected_probs": corrected_all, + "incremental_labels": np.asarray(incremental_labels, dtype=int), + } + + return estimated_prior, details \ No newline at end of file diff --git a/quapy/method/aggregative.py b/quapy/method/aggregative.py index f678546..23341b6 100644 --- a/quapy/method/aggregative.py +++ b/quapy/method/aggregative.py @@ -368,6 +368,65 @@ class BinaryAggregativeQuantifier(AggregativeQuantifier, BinaryQuantifier): return super().fit(X, y) +class ImportanceWeightQuantifier(AggregativeQuantifier, ABC): + """ + Abstract mixin for aggregative quantifiers that estimate the target prevalence by first computing a vector + of importance weights :math:`w_y=Q(y)/P(y)` (with :math:`P` and :math:`Q` the training and target + distributions), and then rescaling the training prevalence by these weights, i.e., + :math:`\\hat{p}(y) \\propto w_y \\cdot P(y)`. + + Subclasses must set a fitted attribute `self.train_prevalence_` (typically in :meth:`aggregation_fit`) and + must implement :meth:`_weights_from_predictions`. This class provides a template implementation of + :meth:`aggregate`, together with :meth:`get_importance_weights` and :meth:`quantify_and_weigh`, none of + which mutate any internal state, so that they are all safe to call concurrently (e.g., from different + threads) on the same fitted instance for different batches of target instances. + """ + + @abstractmethod + def _weights_from_predictions(self, classif_predictions) -> np.ndarray: + """ + Computes the vector of importance weights from the classifier predictions on a batch of (unlabelled) + target instances. + + :param classif_predictions: array-like with the classifier predictions (crisp or soft, depending on + the subclass) for the target instances + :return: np.ndarray of shape `(n_classes,)` + """ + ... + + def _prevalence_from_weights(self, weights: np.ndarray) -> np.ndarray: + estimate = self.train_prevalence_ * weights + return F.normalize_prevalence(estimate, method=self.norm) + + def aggregate(self, classif_predictions): + weights = self._weights_from_predictions(classif_predictions) + return self._prevalence_from_weights(weights) + + def get_importance_weights(self, instances) -> np.ndarray: + """ + Estimates the vector of importance weights :math:`w_y=Q(y)/P(y)` for the given (unlabelled) target + instances. + + :param instances: array-like of shape `(n_instances, n_dimensions)`, the target instances + :return: np.ndarray of shape `(n_classes,)` + """ + classif_predictions = self.classify(instances) + return self._weights_from_predictions(classif_predictions) + + def quantify_and_weigh(self, instances): + """ + Jointly returns the estimated target prevalence and the importance weights used to obtain it, computed + from a single pass of classifier predictions over the given (unlabelled) target instances. + + :param instances: array-like of shape `(n_instances, n_dimensions)`, the target instances + :return: a tuple `(prevalence, weights)`, both np.ndarray of shape `(n_classes,)` + """ + classif_predictions = self.classify(instances) + weights = self._weights_from_predictions(classif_predictions) + prevalence = self._prevalence_from_weights(weights) + return prevalence, weights + + # Methods # ------------------------------------ class CC(AggregativeCrispQuantifier): @@ -456,8 +515,10 @@ class ACC(AggregativeCrispQuantifier): :param str method: adjustment method to be used: - * 'inversion': matrix inversion method based on the matrix equality :math:`P(C)=P(C|Y)P(Y)`, - which tries to invert :math:`P(C|Y)` matrix. + * 'inversion': matrix inversion method. Based on the matrix equality :math:`q=M p`, with + :math:`q` the prevalence vector estimated by CC, :math:`M` the matrix with entries :math:`i,j` representing + :math:`P(\\hat{Y}=i|Y=j)`, and :math:`p` the sought class prevalence vector, the matrix inversion + tries to solve for :math:`p=M^{-1} q` * 'invariant-ratio': invariant ratio estimator of `Vaz et al. 2018 `_, which replaces the last equation with the normalization condition. @@ -467,8 +528,7 @@ class ACC(AggregativeCrispQuantifier): strictly less than `n_classes`. * 'exact-cc': if the matrix is not of full rank, returns `p_c` as the estimates, which corresponds to no adjustment (i.e., the classify and count method. See :class:`quapy.method.aggregative.CC`) - * 'exact': deprecated, defaults to 'exact-cc' - * 'minimize': minimizes the L2 norm of :math:`|Ax-B|`. This one generally works better, and is the + * 'minimize': minimizes the squared L2 norm of :math:`|Ax-B|`. This one generally works better, and is the default parameter. More details about this can be consulted in `Bunse, M. "On Multi-Class Extensions of Adjusted Classify and Count", on proceedings of the 2nd International Workshop on Learning to Quantify: Methods and Applications (LQ 2022), ECML/PKDD 2022, Grenoble (France) @@ -489,7 +549,7 @@ class ACC(AggregativeCrispQuantifier): classifier: BaseEstimator = None, fit_classifier = True, val_split = 5, - solver: Literal['minimize', 'exact', 'exact-raise', 'exact-cc'] = 'minimize', + solver: Literal['minimize', 'exact-raise', 'exact-cc'] = 'minimize', method: Literal['inversion', 'invariant-ratio'] = 'inversion', norm: Literal['clip', 'mapsimplex', 'condsoftmax'] = 'clip', n_jobs=None, @@ -500,7 +560,7 @@ class ACC(AggregativeCrispQuantifier): self.method = method self.norm = norm - SOLVERS = ['exact', 'minimize', 'exact-raise', 'exact-cc'] + SOLVERS = ['minimize', 'exact-raise', 'exact-cc'] METHODS = ['inversion', 'invariant-ratio'] NORMALIZATIONS = ['clip', 'mapsimplex', 'condsoftmax', None] @@ -551,8 +611,8 @@ class ACC(AggregativeCrispQuantifier): @classmethod def getPteCondEstim(cls, classes, y, y_): """ - Estimate the matrix with entry (i,j) being the estimate of P(hat_yi|yj), that is, the probability that a - document that belongs to yj ends up being classified as belonging to yi + Estimate the matrix with entry (i,j) being the estimate of P(hat_yi|yj), that is, the probability that an + instance that belongs to class j ends up being classified as belonging to class i :param classes: array-like with the class names :param y: array-like with the true labels @@ -570,15 +630,234 @@ class ACC(AggregativeCrispQuantifier): return conf def aggregate(self, classif_predictions): - prevs_estim = self.cc.aggregate(classif_predictions) - estimate = F.solve_adjustment( + prevs_estim_cc = self.cc.aggregate(classif_predictions) + estimate = ACC.solve_adjustment( class_conditional_rates=self.Pte_cond_estim_, - unadjusted_counts=prevs_estim, + unadjusted_counts=prevs_estim_cc, solver=self.solver, method=self.method, ) return F.normalize_prevalence(estimate, method=self.norm) + @classmethod + def solve_adjustment(cls, + class_conditional_rates: np.ndarray, + unadjusted_counts: np.ndarray, + method: Literal["inversion", "invariant-ratio"], + solver: Literal["minimize", "exact-raise", "exact-cc"]) -> np.ndarray: + """ + Function that tries to solve for :math:`p` the equation :math:`q = M p`, where :math:`q` is the vector of + `unadjusted counts` (as estimated, e.g., via classify and count) with :math:`q_i` an estimate of + :math:`P(\\hat{Y}=y_i)`, and where :math:`M` is the matrix of `class-conditional rates` with :math:`M_{ij}` an + estimate of :math:`P(\\hat{Y}=y_i|Y=y_j)`. + + :param class_conditional_rates: array of shape `(n_classes, n_classes,)` with entry `(i,j)` being the estimate + of :math:`P(\\hat{Y}=y_i|Y=y_j)`, that is, the probability that an instance ends up being classified as + belonging to class :math:`y_i` given it actually belonged to class :math:`y_j` + + :param unadjusted_counts: array of shape `(n_classes,)` containing the unadjusted prevalence values (e.g., as + estimated by CC or PCC) + + :param str method: indicates the adjustment method to be used. Valid options are: + + * `inversion`: tries to solve the equation :math:`q = M p` as :math:`p = M^{-1} q` where + :math:`M^{-1}` is the matrix inversion of :math:`M`. This inversion may not exist in + degenerated cases. + * `invariant-ratio`: invariant ratio estimator of `Vaz et al. 2018 `_, + which replaces the last equation in :math:`M` with the normalization condition (i.e., that the sum of + all prevalence values must equal 1). + + :param str solver: the method to use for solving the system of linear equations. Valid options are: + + * `exact-raise`: tries to solve the system using matrix inversion. Raises an error if the matrix has rank + strictly lower than `n_classes`. + * `exact-cc`: if the matrix is not full rank, returns :math:`q` (i.e., the unadjusted counts) as the estimates + * `minimize`: minimizes a loss, so the solution always exists + """ + + A = np.asarray(class_conditional_rates, dtype=float).copy() + B = np.asarray(unadjusted_counts, dtype=float).copy() + + if method == "inversion": + pass # leave A and B unchanged + elif method == "invariant-ratio": + # Change the last equation to replace it with the normalization condition; + # copy first so this does not mutate the caller's arrays (np.asarray above + # returns the same object, not a copy, when the input is already float64) + A[-1, :] = 1.0 + B[-1] = 1.0 + else: + raise ValueError(f"unknown {method=}") + + if solver in ["exact-raise", "exact-cc"]: + try: + return np.linalg.solve(A, B) + except np.linalg.LinAlgError: + if solver=='exact-cc': + return unadjusted_counts + else: + raise + elif solver == "minimize": + def loss(prev): + return np.linalg.norm(A @ prev - B) + return F.optim_minimize(loss, n_classes=A.shape[0], return_loss=False) + else: + raise ValueError(f"unknown {solver=}") + + +class BBSEhard(ImportanceWeightQuantifier, AggregativeCrispQuantifier): + """ + `Black Box Shift Estimator` (BBSE) hard aims at finding the importance weights :math:`w_i=Q(i)/P(i)`, + with :math:`P` and :math:`Q` the training and test distributions. + BBSE is similar in spirit to ACC, but it solves the problem :math:`q=C w`, with :math:`q` the prevalence vector + estimated by CC, :math:`C` the matrix with entries :math:`i,j` representing :math:`P(\\hat{Y}=i,Y=j)`, and + :math:`w` the sought vector of importance weights. The `hard` variant estimates these quantities using + crisp counts. + + BBSE was proposed in + `Lipton, Z., Wang, Y. X., & Smola, A. (2018, July). + Detecting and correcting for label shift with black box predictors. + In International conference on machine learning (pp. 3122-3130). PMLR. + `_. + + :param classifier: a scikit-learn's BaseEstimator, or None, in which case the classifier is taken to be + the one indicated in `qp.environ['DEFAULT_CLS']` + + :param fit_classifier: whether to train the learner (default is True). Set to False if the + learner has been trained outside the quantifier. + + :param val_split: specifies the data used for generating classifier predictions. This specification + can be made as float in (0, 1) indicating the proportion of stratified held-out validation set to + be extracted from the training set; or as an integer (default 5), indicating that the predictions + are to be generated in a `k`-fold cross-validation manner (with this integer indicating the value + for `k`); or as a tuple (X,y) defining the specific set of data to use for validation. + + :param str solver: indicates the method to use for solving the system of linear equations. Valid options are: + + * 'exact-raise': tries to solve the system using matrix inversion. Raises an error if the matrix has rank + strictly less than `n_classes`. + * 'exact-cc': if the matrix is not of full rank, returns `p_c` as the estimates, which corresponds to + no adjustment (i.e., the classify and count method. See :class:`quapy.method.aggregative.CC`) + * 'exact': deprecated, defaults to 'exact-cc' + * 'minimize': minimizes the squared L2 norm of :math:`|Ax-B|`. This one generally works better, and is the + default parameter. More details about this can be consulted in + `Tachet des Combes, R., Zhao, H., Wang, Y. X., & Gordon, G. J. (2020). + Domain adaptation with conditional distribution matching and generalized label shift. + Advances in Neural Information Processing Systems, 33, 19276-19289. + `_. + + :param str norm: the method to use for normalization. + + * `clip`, the values are clipped to the range [0,1] and then L1-normalized. + * `mapsimplex` projects vectors onto the probability simplex. This implementation relies on + `Mathieu Blondel's projection_simplex_sort `_ + * `condsoftmax`, applies a softmax normalization only to prevalence vectors that lie outside the simplex + + :param n_jobs: number of parallel workers + """ + + def __init__( + self, + classifier: BaseEstimator = None, + fit_classifier = True, + val_split = 5, + solver: Literal['minimize', 'exact-raise', 'exact-cc'] = 'minimize', + norm: Literal['clip', 'mapsimplex', 'condsoftmax'] = 'clip', + n_jobs=None, + ): + super().__init__(classifier, fit_classifier, val_split) + self.n_jobs = qp._get_njobs(n_jobs) + self.solver = solver + self.norm = norm + + def _check_init_parameters(self): + if self.solver not in ACC.SOLVERS: + raise ValueError(f"unknown solver; valid ones are {ACC.SOLVERS}") + if self.norm not in ACC.NORMALIZATIONS: + raise ValueError(f"unknown normalization; valid ones are {ACC.NORMALIZATIONS}") + + def aggregation_fit(self, classif_predictions, labels): + """ + Estimates the misclassification rates. + :param classif_predictions: array-like with the predicted labels + :param labels: array-like with the true labels associated to each predicted label + """ + true_labels = labels + pred_labels = classif_predictions + self.cc = CC(self.classifier, fit_classifier=False) + self.confusion = BBSEhard.getConfusionJointProb(self.classifier.classes_, true_labels, pred_labels) + self.train_prevalence_ = F.prevalence_from_labels(labels, classes=self.classifier.classes_) + + @classmethod + def getConfusionJointProb(cls, classes, y, y_): + """ + Estimate the matrix with entry (i,j) being the estimate of :math:`P(\\hat{Y}=i,Y=j)` + + :param classes: array-like with the class names + :param y: array-like with the true labels + :param y_: array-like with the estimated labels + :return: np.ndarray + """ + conf = confusion_matrix(y, y_, labels=classes).T + joint_probs = conf / conf.sum() + return joint_probs + + def _weights_from_predictions(self, classif_predictions): + prevs_estim_cc = self.cc.aggregate(classif_predictions) + weights = BBSEhard.solve_importance_weights( + joint_probs=self.confusion, + unadjusted_counts=prevs_estim_cc, + solver=self.solver, + ) + return np.clip(weights, 0.0, None) + + @classmethod + def solve_importance_weights(cls, + joint_probs: np.ndarray, + unadjusted_counts: np.ndarray, + solver: Literal["minimize", "exact-raise", "exact-cc"]) -> np.ndarray: + """ + Function that tries to solve for :math:`p` the equation :math:`q = C w`, where :math:`q` is the vector of + `unadjusted counts` (as estimated, e.g., via classify and count) with :math:`q_i` an estimate of + :math:`P(\\hat{Y}=y_i)`, and where :math:`C` is the confusion matrix distribution with :math:`C_{ij}` an + estimate of :math:`P(\\hat{Y}=y_i,Y=y_j)`. + + :param joint_probs: array of shape `(n_classes, n_classes,)` with entry `(i,j)` being the estimate + of :math:`P(\\hat{Y}=y_i,Y=y_j)` + + :param unadjusted_counts: array of shape `(n_classes,)` containing the unadjusted prevalence values (e.g., as + estimated by CC or PCC) + + :param str solver: the method to use for solving the system of linear equations. Valid options are: + + * `exact-raise`: tries to solve the system using matrix inversion. Raises an error if the matrix has rank + strictly lower than `n_classes`. + * `exact-cc`: if the matrix is not full rank, returns a vector of ones as the weights + * `minimize`: minimizes a loss, so the solution always exists + """ + A = np.asarray(joint_probs, dtype=float) + B = np.asarray(unadjusted_counts, dtype=float) + all_ones = np.full_like(B, fill_value=1., dtype=float) + + if solver in ["exact-raise", "exact-cc"]: + try: + return np.linalg.solve(A, B) + except np.linalg.LinAlgError: + if solver=='exact-cc': + return all_ones + else: + raise + elif solver == "minimize": + def loss(prev): + return np.linalg.norm(A @ prev - B) + + n_dims = len(all_ones) + bounds = [(0, np.inf)] * n_dims + return F.optim_minimize(loss, n_classes=A.shape[0], x0=all_ones, bounds=bounds, constraints=()) + else: + raise ValueError(f"unknown {solver=}") + + class PACC(AggregativeSoftQuantifier): """ @@ -615,7 +894,6 @@ class PACC(AggregativeSoftQuantifier): Raises an error if the matrix has rank strictly less than `n_classes`. * 'exact-cc': if the matrix is not of full rank, returns `p_c` as the estimates, which corresponds to no adjustment (i.e., the classify and count method. See :class:`quapy.method.aggregative.CC`) - * 'exact': deprecated, defaults to 'exact-cc' * 'minimize': minimizes the L2 norm of :math:`|Ax-B|`. This one generally works better, and is the default parameter. More details about this can be consulted in `Bunse, M. "On Multi-Class Extensions of Adjusted Classify and Count", on proceedings of the 2nd International Workshop on Learning to @@ -637,7 +915,7 @@ class PACC(AggregativeSoftQuantifier): classifier: BaseEstimator = None, fit_classifier=True, val_split=5, - solver: Literal['minimize', 'exact', 'exact-raise', 'exact-cc'] = 'minimize', + solver: Literal['minimize', 'exact-raise', 'exact-cc'] = 'minimize', method: Literal['inversion', 'invariant-ratio'] = 'inversion', norm: Literal['clip', 'mapsimplex', 'condsoftmax'] = 'clip', n_jobs=None @@ -671,7 +949,7 @@ class PACC(AggregativeSoftQuantifier): def aggregate(self, classif_posteriors): prevs_estim = self.pcc.aggregate(classif_posteriors) - estimate = F.solve_adjustment( + estimate = ACC.solve_adjustment( class_conditional_rates=self.Pte_cond_estim_, unadjusted_counts=prevs_estim, solver=self.solver, @@ -693,7 +971,120 @@ class PACC(AggregativeSoftQuantifier): return confusion.T -class RLLS(AggregativeSoftQuantifier): +class BBSEsoft(ImportanceWeightQuantifier, AggregativeSoftQuantifier): + """ + `Black Box Shift Estimator` (BBSE) soft, the probabilistic variant of :class:`BBSEhard` that relies on the + posterior probabilities returned by a probabilistic classifier, instead of on crisp counts, to estimate the + joint distribution :math:`P(\\hat{Y}=i,Y=j)`. As in :class:`BBSEhard`, the sought importance weights + :math:`w_i=Q(i)/P(i)` (with :math:`P` and :math:`Q` the training and test distributions) are obtained by + solving :math:`q=C w`, with :math:`q` the (now probabilistic) prevalence vector estimated by PCC and + :math:`C` the joint-probability matrix. + + BBSE was proposed in + `Lipton, Z., Wang, Y. X., & Smola, A. (2018, July). + Detecting and correcting for label shift with black box predictors. + In International conference on machine learning (pp. 3122-3130). PMLR. + `_. + + :param classifier: a scikit-learn's BaseEstimator, or None, in which case the classifier is taken to be + the one indicated in `qp.environ['DEFAULT_CLS']` + + :param fit_classifier: whether to train the learner (default is True). Set to False if the + learner has been trained outside the quantifier. + + :param val_split: specifies the data used for generating classifier predictions. This specification + can be made as float in (0, 1) indicating the proportion of stratified held-out validation set to + be extracted from the training set; or as an integer (default 5), indicating that the predictions + are to be generated in a `k`-fold cross-validation manner (with this integer indicating the value + for `k`); or as a tuple (X,y) defining the specific set of data to use for validation. + + :param str solver: indicates the method to use for solving the system of linear equations. Valid options are: + + * 'exact-raise': tries to solve the system using matrix inversion. Raises an error if the matrix has rank + strictly less than `n_classes`. + * 'exact-cc': if the matrix is not of full rank, returns `p_c` as the estimates, which corresponds to + no adjustment (i.e., the classify and count method. See :class:`quapy.method.aggregative.CC`) + * 'minimize': minimizes the squared L2 norm of :math:`|Ax-B|`. This one generally works better, and is the + default parameter. More details about this can be consulted in + `Tachet des Combes, R., Zhao, H., Wang, Y. X., & Gordon, G. J. (2020). + Domain adaptation with conditional distribution matching and generalized label shift. + Advances in Neural Information Processing Systems, 33, 19276-19289. + `_. + + :param str norm: the method to use for normalization. + + * `clip`, the values are clipped to the range [0,1] and then L1-normalized. + * `mapsimplex` projects vectors onto the probability simplex. This implementation relies on + `Mathieu Blondel's projection_simplex_sort `_ + * `condsoftmax`, applies a softmax normalization only to prevalence vectors that lie outside the simplex + + :param n_jobs: number of parallel workers + """ + + def __init__( + self, + classifier: BaseEstimator = None, + fit_classifier = True, + val_split = 5, + solver: Literal['minimize', 'exact-raise', 'exact-cc'] = 'minimize', + norm: Literal['clip', 'mapsimplex', 'condsoftmax'] = 'clip', + n_jobs=None, + ): + super().__init__(classifier, fit_classifier, val_split) + self.n_jobs = qp._get_njobs(n_jobs) + self.solver = solver + self.norm = norm + + def _check_init_parameters(self): + if self.solver not in ACC.SOLVERS: + raise ValueError(f"unknown solver; valid ones are {ACC.SOLVERS}") + if self.norm not in ACC.NORMALIZATIONS: + raise ValueError(f"unknown normalization; valid ones are {ACC.NORMALIZATIONS}") + + def aggregation_fit(self, classif_predictions, labels): + """ + Estimates the joint distribution P(hat_Y,Y), using posterior probabilities in place of crisp counts. + + :param classif_predictions: array-like with posterior probabilities + :param labels: array-like with the true labels associated to each vector of posterior probabilities + """ + posteriors = classif_predictions + true_labels = labels + self.pcc = PCC(self.classifier, fit_classifier=False) + self.confusion = BBSEsoft.getConfusionJointProb(self.classifier.classes_, true_labels, posteriors) + self.train_prevalence_ = F.prevalence_from_labels(labels, classes=self.classifier.classes_) + + @classmethod + def getConfusionJointProb(cls, classes, y, posteriors): + """ + Estimate the matrix with entry (i,j) being the estimate of :math:`P(\\hat{Y}=i,Y=j)`, using the + posterior probabilities of the instances belonging to class :math:`j` in place of their hard predictions. + + :param classes: array-like with the class names + :param y: array-like with the true labels + :param posteriors: array-like of shape `(n_instances, n_classes,)` with posterior probabilities + :return: np.ndarray + """ + n_classes = len(classes) + joint_probs = np.zeros((n_classes, n_classes), dtype=float) + for j, class_ in enumerate(classes): + idx = y == class_ + if idx.any(): + joint_probs[:, j] = posteriors[idx].sum(axis=0) + joint_probs /= joint_probs.sum() + return joint_probs + + def _weights_from_predictions(self, classif_posteriors): + prevs_estim_pcc = self.pcc.aggregate(classif_posteriors) + weights = BBSEhard.solve_importance_weights( + joint_probs=self.confusion, + unadjusted_counts=prevs_estim_pcc, + solver=self.solver, + ) + return np.clip(weights, 0.0, None) + + +class RLLS(ImportanceWeightQuantifier, AggregativeSoftQuantifier): """ `Regularized Learning for Domain Adaptation under Label Shifts `_, used here as an aggregative @@ -750,7 +1141,6 @@ class RLLS(AggregativeSoftQuantifier): self.delta = delta self.clip_weights = clip_weights self.norm = norm - self.last_w_ = None def _check_init_parameters(self): _get_cvxpy() @@ -781,18 +1171,15 @@ class RLLS(AggregativeSoftQuantifier): self.pz_ = _rlls_predicted_marginal(classif_predictions, mode=self.mode) self.rho_ = _rlls_compute_3deltaC(len(self.classes_), len(labels), self.delta) - def aggregate(self, classif_posteriors): + def _weights_from_predictions(self, classif_posteriors): qz = _rlls_predicted_marginal(classif_posteriors, mode=self.mode) - w = _rlls_compute_weights( + return _rlls_compute_weights( self.C_zy_, qz, self.pz_, rho=self.alpha * self.rho_, clip=self.clip_weights, ) - self.last_w_ = w - estimate = self.train_prevalence_ * w - return F.normalize_prevalence(estimate, method=self.norm) class EMQ(AggregativeSoftQuantifier): @@ -1831,6 +2218,10 @@ KDEyML = _kdey.KDEyML KDEyHD = _kdey.KDEyHD KDEyCS = _kdey.KDEyCS +from . import _liep + +LEIP = _liep.LEIP + # --------------------------------------------------------------- # aliases diff --git a/quapy/tests/test_methods.py b/quapy/tests/test_methods.py index dfd3b12..7a225dd 100644 --- a/quapy/tests/test_methods.py +++ b/quapy/tests/test_methods.py @@ -8,7 +8,7 @@ from sklearn.linear_model import LogisticRegression from quapy.method import AGGREGATIVE_METHODS, BINARY_METHODS, NON_AGGREGATIVE_METHODS from quapy.method.non_aggregative import DMx, EDx, HDx -from quapy.method.aggregative import ACC, DMy, EDy, KDEyCS, RLLS +from quapy.method.aggregative import ACC, BBSEhard, BBSEsoft, DMy, EDy, KDEyCS, LEIP, RLLS from quapy.method.meta import Ensemble from quapy.functional import check_prevalence_vector from quapy.tests._synthetic import make_dataset @@ -21,6 +21,7 @@ OPTIONAL_AGGREGATIVE_METHODS = { 'PQ', 'RLLS', 'EDy', + 'LEIP', } OPTIONAL_NON_AGGREGATIVE_METHODS = { @@ -178,6 +179,29 @@ class TestMethods(unittest.TestCase): self.assertTrue(check_prevalence_vector(estim_prevalences)) + def test_leip(self): + dataset = TestMethods.tiny_dataset_multiclass + q = LEIP(LogisticRegression(max_iter=2000), val_split=3) + q.fit(*dataset.training.Xy) + estim_prevalences = q.predict(dataset.test.X) + self.assertTrue(check_prevalence_vector(estim_prevalences)) + + def test_leip_fixed_tau(self): + dataset = TestMethods.tiny_dataset_binary + q = LEIP(LogisticRegression(max_iter=2000), val_split=None, tau=0.6) + q.fit(*dataset.training.Xy) + estim_prevalences = q.predict(dataset.test.X) + self.assertTrue(check_prevalence_vector(estim_prevalences)) + + def test_bbse(self): + dataset = TestMethods.tiny_dataset_multiclass + for cls in (BBSEhard, BBSEsoft): + for solver in ('minimize', 'exact-raise', 'exact-cc'): + q = cls(LogisticRegression(max_iter=2000), val_split=3, solver=solver) + q.fit(*dataset.training.Xy) + estim_prevalences = q.predict(dataset.test.X) + self.assertTrue(check_prevalence_vector(estim_prevalences)) + def test_edy(self): try: import quadprog # noqa: F401 From d6fbd13ecd0f17f89237ff73712a5efba8c7bfeb Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Tue, 8 Sep 2026 14:18:41 +0200 Subject: [PATCH 2/9] updated todo --- TODO.txt | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/TODO.txt b/TODO.txt index 12f1f59..fd720cc 100644 --- a/TODO.txt +++ b/TODO.txt @@ -21,4 +21,6 @@ scale each value by per-class thresholds, i.e., [0.33*0.1, 0.33*1, 0.33*1]/sum. check https://docs.python.org/3/library/warnings.html#temporarily-suppressing-warnings - [TODO] nmd and md are not selectable from qp.evaluation.evaluate as a string - [TODO] add https://www.kaggle.com/datasets/Cornell-University/arxiv dataset? in the paper "Online Adaptation to - Label Distribution Shift" they say it provides a natural label shift over time \ No newline at end of file + Label Distribution Shift" they say it provides a natural label shift over time +- [TODO] add Bhattacharyya distance (https://en.wikipedia.org/wiki/Bhattacharyya_distance) + (which is actually not a proper distance) \ No newline at end of file From 89548d3a8bef7dd4d0e62011f39e672bbeac64b4 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Tue, 8 Sep 2026 14:47:03 +0200 Subject: [PATCH 3/9] Add GMNet, a Gaussian-mixture neural quantifier Ports GMNet (from https://github.com/pglez84/gmnet) into quapy/method/_gmnet.py, mirroring how HistNetQ was ported: dropping that repo's quantificationlib-backed bag generators in favor of QuaPy's own sampling protocols, and adding geotorch (now a 'neural' extra dependency) to keep the Gaussian layers' covariance matrices positive-definite during training. - GMNet represents each bag instance by its likelihood under one or more learned mixtures of Gaussians ("GM branches"), mean-pools these representations over the bag, and predicts prevalence from the result. Supports multiple stacked GM branches with an optional CKA-regularization term encouraging their latent representations to be dissimilar. - Fixes two aspects of the original architecture that assumed a fixed, training-time bag_size baked into the network (a reshape step, and forward-hook-based activation capture for CKA): both are now computed from the actual input shape/plain attributes at forward time, so the model also works on predict()'s arbitrary-sized test samples, not just same-size bags. - Factors the bag-based training loop shared by HistNetQ and GMNet (bag generation, fit/fit_from_samples, early stopping, LR scheduling, checkpointing, predict) out of _histnet.py into a new BagTrainedQuantifier base class in quapy/method/_neural_bags.py; HistNetQ's public API and behavior are unchanged. - Aliased in meta.py (torch/geotorch-optional, mirroring HistNetQ/QuaNet) and registered in META_METHODS. - Adds test_gmnet covering single-branch and multi-branch+CKA (via fit_from_samples/mix_bags) variants. Co-Authored-By: Claude Sonnet 5 --- quapy/method/__init__.py | 3 +- quapy/method/_gmnet.py | 328 ++++++++++++++++++++++++++ quapy/method/_histnet.py | 355 ++++------------------------ quapy/method/_neural_bags.py | 443 +++++++++++++++++++++++++++++++++++ quapy/method/meta.py | 11 + quapy/tests/test_methods.py | 32 +++ setup.py | 2 +- 7 files changed, 864 insertions(+), 310 deletions(-) create mode 100644 quapy/method/_gmnet.py create mode 100644 quapy/method/_neural_bags.py diff --git a/quapy/method/__init__.py b/quapy/method/__init__.py index 695ac5f..d295382 100644 --- a/quapy/method/__init__.py +++ b/quapy/method/__init__.py @@ -81,7 +81,8 @@ NON_AGGREGATIVE_METHODS = { META_METHODS = { meta.Ensemble, meta.QuaNet, - meta.HistNetQ + meta.HistNetQ, + meta.GMNet } QUANTIFICATION_METHODS = AGGREGATIVE_METHODS | NON_AGGREGATIVE_METHODS | META_METHODS diff --git a/quapy/method/_gmnet.py b/quapy/method/_gmnet.py new file mode 100644 index 0000000..195681d --- /dev/null +++ b/quapy/method/_gmnet.py @@ -0,0 +1,328 @@ +""" +GMNet implementation. + +Ported from the reference implementation at https://github.com/pglez84/gmnet (the `GMNet`/ +`DLQuantification` classes in that repo), adapted to QuaPy's own protocol-based sample generation +(replacing that repo's custom, `quantificationlib`-backed bag generators), and reusing the shared +bag-based training loop already factored out for :class:`quapy.method._histnet.HistNetQ` (see +:class:`quapy.method._neural_bags.BagTrainedQuantifier`). + +The overall architecture is: one or more "GM branches" -- each a small per-branch feature extractor +followed by a layer of Gaussian likelihoods (a :class:`_GMLayer`) evaluated at every instance of a bag +-- concatenated and mean-pooled over the bag, followed by the shared quantification MLP head. Like +HistNetQ (and QuaNet), GMNet is trained end-to-end by minimizing a quantification loss over samples +("bags") of known prevalence, rather than over individually labeled instances. + +Two deliberate deviations from the reference implementation, both required for the model to satisfy +QuaPy's `predict(X)` contract (i.e., to be usable on a real test collection of arbitrary size, as +opposed to only on bags resampled at the fixed `bag_size` used for training): + +* the original `GMNet_Module` reshapes each branch's per-instance likelihoods around a *fixed*, + constructor-time `bag_size` (via `torch.nn.Unflatten(0, (-1, bag_size))`), which only works when + every forward pass is fed bags of exactly that size. Here, the reshape is instead computed from the + actual input shape at forward time (see :class:`_GMBranch`), which is equivalent when the bag size + matches but also supports bags (or, at prediction time, whole test samples) of any other size. +* the forward hooks used by the original code to capture each branch's pre-Gaussian latent activations + (for the CKA regularization term) are replaced by simply storing that activation as an attribute + during `forward` (see :attr:`_GMBranch.latent_activation`), since branches are now implemented with a + plain `forward` method rather than an opaque `torch.nn.Sequential`. +""" +import numpy as np +import scipy.spatial.distance +import torch +import torch.nn as nn + +import geotorch + +from quapy.method._neural_bags import BagTrainedQuantifier +from quapy.protocol import UPP + + +def _cka(latent_activations): + """Feature-space linear CKA (Centered Kernel Alignment), averaged over every pair of latent + activations, following the `CKARegularization` class in the reference implementation. Used to + encourage the Gaussian components learned by different GM branches to capture complementary + (dissimilar) aspects of the instances. + + :param latent_activations: a list of tensors, one per GM branch, all of shape (n_instances, dim_i) + (dim_i may differ across branches) + """ + cka_sum = 0. + n_pairs = 0 + for i in range(len(latent_activations)): + for j in range(i + 1, len(latent_activations)): + x = latent_activations[i] + y = latent_activations[j] + x = x - torch.mean(x, dim=0, keepdim=True) + y = y - torch.mean(y, dim=0, keepdim=True) + dot_product_similarity = torch.norm(torch.matmul(x.t(), y)) ** 2 + normalization_x = torch.norm(torch.matmul(x.t(), x)) + normalization_y = torch.norm(torch.matmul(y.t(), y)) + cka_sum = cka_sum + dot_product_similarity / (normalization_x * normalization_y) + n_pairs += 1 + return cka_sum / n_pairs + + +class _GMLayer(nn.Module): + """A layer of `num_gaussians` (unnormalized) Gaussian likelihoods, evaluated at every instance of + a bag. `centers` and `covariance` are learned; `covariance` is constrained to stay positive-definite + throughout training via `geotorch.positive_definite`. + """ + + def __init__(self, n_features, num_gaussians): + super().__init__() + self.n_features = n_features + self.num_gaussians = num_gaussians + self.centers = nn.Parameter(torch.rand(num_gaussians, n_features)) + self.covariance = nn.Parameter(torch.eye(n_features).repeat(num_gaussians, 1, 1)) + geotorch.positive_definite(self, "covariance") + + # initialize the centers' covariance from the (squared, halved) nearest-neighbor distance + # between the randomly initialized centers, so that gaussians start with a sensible spread + centers = self.centers.detach().cpu().numpy() + distances = scipy.spatial.distance.cdist(centers, centers) + np.fill_diagonal(distances, np.inf) + cov = (np.mean(np.min(distances, axis=1)) / 2) ** 2 + self.covariance = torch.eye(n_features).repeat(num_gaussians, 1, 1) * cov + + def forward(self, x): + # x: (batch_size, bag_size, n_features) + centers = self.centers.unsqueeze(0).unsqueeze(0) # (1, 1, num_gaussians, n_features) + diff = x.unsqueeze(2) - centers # (batch_size, bag_size, num_gaussians, n_features) + + cov_inv = torch.inverse(self.covariance) + det_cov = torch.linalg.det(self.covariance) + + mahalanobis = torch.einsum('...i,...ij,...j->...', diff, cov_inv.unsqueeze(0).unsqueeze(0), diff) + normalization_term = torch.log((2 * torch.pi) ** self.n_features * det_cov).unsqueeze(0).unsqueeze(0) + log_probs = -0.5 * (mahalanobis + normalization_term) + return torch.exp(log_probs) # (batch_size, bag_size, num_gaussians) + + +class _GMBranch(nn.Module): + """One GM branch: an optional small MLP mapping the (already feature-extracted) instances into a + `gaussian_dimensions`-sized latent space, followed by a Sigmoid, a :class:`_GMLayer`, and a + BatchNorm applied instance-wise (i.e., over the merged batch*bag_size dimension, matching the + reference implementation). + """ + + def __init__(self, input_size, num_gaussians, gaussian_dimensions, hidden_size_fe, dropout_fe): + super().__init__() + self.pre = nn.Sequential() + prev_size = input_size + latent_size = gaussian_dimensions if gaussian_dimensions is not None else input_size + if gaussian_dimensions is not None: + for j, layer_size in enumerate(hidden_size_fe or ()): + self.pre.add_module(f'hidden_{j}', nn.Linear(prev_size, layer_size)) + self.pre.add_module(f'leakyrelu_{j}', nn.LeakyReLU()) + self.pre.add_module(f'dropout_{j}', nn.Dropout(dropout_fe)) + prev_size = layer_size + self.pre.add_module('latent_linear', nn.Linear(prev_size, gaussian_dimensions)) + self.pre.add_module('sigmoid', nn.Sigmoid()) + + self.gm_layer = _GMLayer(n_features=latent_size, num_gaussians=num_gaussians) + self.batch_norm = nn.BatchNorm1d(num_features=num_gaussians) + self.output_size = num_gaussians + self.latent_activation = None # populated on every forward(), read by GMNet's CKA regularization + + def forward(self, x): + # x: (batch_size, bag_size, input_size) + batch_size, bag_size = x.shape[0], x.shape[1] + latent = self.pre(x) + self.latent_activation = latent.reshape(-1, latent.shape[-1]) + likelihoods = self.gm_layer(latent) # (batch_size, bag_size, num_gaussians) + flat = self.batch_norm(likelihoods.reshape(batch_size * bag_size, -1)) + return flat.reshape(batch_size, bag_size, -1) + + +class _GMNetModule(nn.Module): + """The quantification module for GMNet: one or more :class:`_GMBranch` instances, each producing a + per-instance representation that is concatenated across branches and mean-pooled over the bag, as + required by :class:`quapy.method._neural_bags.BagTrainedQuantifier`. + """ + + def __init__(self, input_size, num_gaussians, n_gm_layers, gaussian_dimensions, hidden_size_fe=None, + dropout_fe=0., cka_regularization=0.): + super().__init__() + if len(num_gaussians) != n_gm_layers: + raise ValueError('num_gaussians should be a tuple of the same size as n_gm_layers') + if len(gaussian_dimensions) != n_gm_layers: + raise ValueError('gaussian_dimensions should be a tuple of the same size as n_gm_layers') + + self.n_gm_layers = n_gm_layers + self.cka_regularization = cka_regularization + self.branches = nn.ModuleList([ + _GMBranch(input_size, num_gaussians[i], gaussian_dimensions[i], hidden_size_fe, dropout_fe) + for i in range(n_gm_layers) + ]) + self.output_size = sum(num_gaussians) + + def forward(self, x): + outputs = [branch(x) for branch in self.branches] + return torch.mean(torch.cat(outputs, dim=-1), dim=1) + + def apply_regularization(self): + """Whether the CKA regularization term should be added to the training loss: requires at least + two GM branches (CKA is a pairwise measure) and a nonzero `cka_regularization` weight.""" + return self.n_gm_layers > 1 and self.cka_regularization != 0 + + def regularization_term(self): + latent_activations = [branch.latent_activation for branch in self.branches] + return self.cka_regularization * _cka(latent_activations) + + +class GMNet(BagTrainedQuantifier): + """ + Implementation of `GMNet `_, a neural network for quantification + that represents each instance of a bag by its likelihood under one or more learned mixtures of + Gaussians, mean-pools these representations over the bag, and predicts the class prevalence from the + result, trained end-to-end by minimizing a quantification loss over many samples ("bags") of known + prevalence. + + Like :class:`quapy.method._histnet.HistNetQ` and :class:`quapy.method.meta.QuaNet`, GMNet does not + follow the classify-then-aggregate pattern of :class:`quapy.method.aggregative.AggregativeQuantifier`; + it is instead trained and evaluated end-to-end on whole bags (see + :class:`quapy.method._neural_bags.BagTrainedQuantifier` for the shared training/prediction logic, + including the two entry points, :meth:`fit` and :meth:`fit_from_samples`). + + :param feature_extraction_module: a `torch.nn.Module` exposing an `output_size` attribute, used to + embed each instance before it is passed to every GM branch. If None (default), an identity + module is used, i.e., the instances in `X` are assumed to already be in their final numeric + representation. + :param n_gm_layers: number of GM branches (default 1). + :param num_gaussians: number of gaussians per branch: either a single int (used for every branch) or + a tuple/list of `n_gm_layers` ints (default 4). + :param gaussian_dimensions: dimensionality of the latent space in which each branch's gaussians live: + either a single int/None (used for every branch) or a tuple/list of `n_gm_layers` int/None + values. If None for a given branch, that branch's gaussians operate directly on the + feature-extracted instances, with no extra per-branch projection (default None). + :param hidden_size_fe: sizes of the hidden layers of the small per-branch MLP that maps the + feature-extracted instances into the latent space (only used when `gaussian_dimensions` is not + None for the corresponding branch); default None (no hidden layers, i.e., a single linear + projection). + :param dropout_fe: dropout applied after each of the `hidden_size_fe` layers (default 0). + :param cka_regularization: weight of the CKA regularization term encouraging the different branches' + latent representations to be dissimilar; only applied when `n_gm_layers > 1` (default 0, i.e., + disabled). + :param linear_sizes: tuple of ints with the sizes of the linear layers used in the shared + quantification head, after the GM branches (default empty, i.e., only the final classification + layer is used). + :param dropout: dropout applied after each of the `linear_sizes` layers (default 0). + :param output_function: either 'softmax' or 'normalize' (L1); both yield a valid prevalence vector + (default 'softmax'). + :param bag_size: number of instances per training/validation bag (default 500). + :param n_bags_train: number of bags generated per training epoch (default 500). + :param n_bags_val: number of bags generated per validation epoch (default 500). + :param train_epochs: maximum number of training epochs (default 200). + :param patience: number of epochs without improvement in validation loss before early-stopping + (default 20). + :param start_lr: initial learning rate (default 1e-3). + :param end_lr: once the learning rate decays below this value, training stops (default 1e-6). + :param lr_factor: factor by which the learning rate is reduced after `patience` epochs without + improvement (default 0.1). + :param weight_decay: L2 regularization (default 0). + :param quant_loss: the quantification loss to minimize (default `torch.nn.L1Loss()`), called as + `quant_loss(true_prevalences, predicted_prevalences)`. + :param batch_size: number of bags per gradient update (default 16). + :param protocol: the :class:`quapy.protocol.AbstractStochasticSeededProtocol` subclass used by + :meth:`fit` to resample bags from the given labelled collection (default + :class:`quapy.protocol.UPP`, which draws bags with prevalence sampled uniformly at random from + the simplex). + :param protocol_params: dict of extra keyword arguments passed to `protocol` (besides `data`, + `sample_size`, `repeats`, and `random_state`, which are set internally); default None. + :param val_split: float in (0,1), the proportion of the collection given to :meth:`fit` that is held + out (via stratified sampling) for validation and early stopping (default 0.4). + :param device: `'cpu'` or `'cuda'` (default 'cpu'). + :param random_state: seed used for the train/validation split and for the (fixed) validation + sampling sequence, as well as for the random initialization of the GM branches (default 0). + :param checkpointdir: directory where the best model found during training is stored (default + '../checkpoint'). + :param checkpointname: name of the checkpoint file; if None (default), a random name is generated. + :param verbose: verbosity level; if >0, shows a progress bar with the current losses (default 0). + """ + + def __init__(self, + feature_extraction_module=None, + n_gm_layers=1, + num_gaussians=4, + gaussian_dimensions=None, + hidden_size_fe=None, + dropout_fe=0., + cka_regularization=0., + linear_sizes=(), + dropout=0., + output_function='softmax', + bag_size=500, + n_bags_train=500, + n_bags_val=500, + train_epochs=200, + patience=20, + start_lr=1e-3, + end_lr=1e-6, + lr_factor=0.1, + weight_decay=0., + quant_loss=None, + batch_size=16, + protocol=UPP, + protocol_params=None, + val_split=0.4, + device='cpu', + random_state=0, + checkpointdir='../checkpoint', + checkpointname=None, + verbose=0): + super().__init__( + feature_extraction_module=feature_extraction_module, + linear_sizes=linear_sizes, + dropout=dropout, + output_function=output_function, + bag_size=bag_size, + n_bags_train=n_bags_train, + n_bags_val=n_bags_val, + train_epochs=train_epochs, + patience=patience, + start_lr=start_lr, + end_lr=end_lr, + lr_factor=lr_factor, + weight_decay=weight_decay, + quant_loss=quant_loss, + batch_size=batch_size, + protocol=protocol, + protocol_params=protocol_params, + val_split=val_split, + device=device, + random_state=random_state, + checkpointdir=checkpointdir, + checkpointname=checkpointname, + verbose=verbose, + ) + self.n_gm_layers = n_gm_layers + self.num_gaussians = num_gaussians if isinstance(num_gaussians, (tuple, list)) \ + else [num_gaussians] * n_gm_layers + self.gaussian_dimensions = gaussian_dimensions if isinstance(gaussian_dimensions, (tuple, list)) \ + else [gaussian_dimensions] * n_gm_layers + self.hidden_size_fe = hidden_size_fe + self.dropout_fe = dropout_fe + self.cka_regularization = cka_regularization + + @property + def _checkpoint_prefix(self): + return 'GMNet' + + def _build_quantmodule(self, n_features): + torch.manual_seed(self.random_state) + return _GMNetModule( + input_size=n_features, + num_gaussians=self.num_gaussians, + n_gm_layers=self.n_gm_layers, + gaussian_dimensions=self.gaussian_dimensions, + hidden_size_fe=self.hidden_size_fe, + dropout_fe=self.dropout_fe, + cka_regularization=self.cka_regularization, + ) + + def _extra_loss(self): + quantmodule = self.model.quantmodule + if quantmodule.apply_regularization(): + return quantmodule.regularization_term() + return 0. diff --git a/quapy/method/_histnet.py b/quapy/method/_histnet.py index ac4cc40..ad58725 100644 --- a/quapy/method/_histnet.py +++ b/quapy/method/_histnet.py @@ -15,33 +15,15 @@ arXiv preprint arXiv:2012.06311 (2020). The overall architecture is: feature_extraction -> Sigmoid -> histogram layer -> small MLP -> softmax, trained by minimizing a quantification loss over samples ("bags") of known prevalence, rather than -over individually labeled instances (in the spirit of QuaNet, see method/_quanet.py). +over individually labeled instances (in the spirit of QuaNet, see method/_quanet.py). The bag-based +training loop itself (bag generation, early stopping, LR scheduling, checkpointing, prediction) is +shared with :class:`quapy.method._gmnet.GMNet` via :class:`quapy.method._neural_bags.BagTrainedQuantifier`. """ -import copy -import os -import random - -import numpy as np import torch import torch.nn as nn -from tqdm import tqdm -from quapy.data import LabelledCollection -from quapy.method.base import BaseQuantifier -from quapy.protocol import AbstractProtocol, UPP -from quapy.util import EarlyStop - - -class _IdentityFeatureExtractionModule(nn.Module): - """Used when no feature extraction module is provided: instances are assumed to already be in - their final numeric representation.""" - - def __init__(self, input_size): - super().__init__() - self.output_size = input_size - - def forward(self, x): - return x +from quapy.method._neural_bags import BagTrainedQuantifier +from quapy.protocol import UPP class _HardHistogramLayer(nn.Module): @@ -99,84 +81,22 @@ class _HardHistogramLayer(nn.Module): return result -class _HistNetModule(nn.Module): - """The full HistNetQ network: feature extraction, histogram, and the quantification MLP.""" +class _SigmoidHistogram(nn.Module): + """The quantification module for HistNetQ: squashes the (already feature-extracted) instances + through a Sigmoid and builds a differentiable histogram of them, as required by + :class:`quapy.method._neural_bags.BagTrainedQuantifier`.""" - def __init__(self, feature_extraction_module, n_classes, n_bins=8, quantiles=False, linear_sizes=(), - dropout=0., output_function='softmax'): + def __init__(self, n_features, n_bins=8, quantiles=False): super().__init__() - self.feature_extraction_module = feature_extraction_module self.sigmoid = nn.Sigmoid() - self.histogram = _HardHistogramLayer( - n_features=feature_extraction_module.output_size, n_bins=n_bins, quantiles=quantiles - ) + self.histogram = _HardHistogramLayer(n_features=n_features, n_bins=n_bins, quantiles=quantiles) + self.output_size = self.histogram.output_size - self.output_function = output_function - self.output_module = nn.Sequential() - prev_size = self.histogram.output_size - for i, linear_size in enumerate(linear_sizes): - self.output_module.add_module(f'linear_{i}', nn.Linear(prev_size, linear_size)) - self.output_module.add_module(f'leakyrelu_{i}', nn.LeakyReLU()) - self.output_module.add_module(f'dropout_{i}', nn.Dropout(dropout)) - prev_size = linear_size - self.output_module.add_module('last_linear', nn.Linear(prev_size, n_classes)) - if output_function == 'softmax': - self.output_module.add_module('softmax', nn.Softmax(dim=1)) - elif output_function == 'normalize': - self.output_module.add_module('relu', nn.ReLU()) - else: - raise ValueError(f"unknown {output_function=}; valid ones are 'softmax', 'normalize'") - - def forward(self, bag): - # bag: (batch_size, bag_size, n_features) - features = self.feature_extraction_module(bag) - features = self.sigmoid(features) - histogram = self.histogram(features) - out = self.output_module(histogram) - if self.output_function == 'normalize': - out = nn.functional.normalize(out, p=1, dim=1) - return out + def forward(self, input): + return self.histogram(self.sigmoid(input)) -def _to_tensor(x, device): - if torch.is_tensor(x): - return x.to(device=device, dtype=torch.float32) - if hasattr(x, 'toarray'): # scipy sparse - x = x.toarray() - return torch.as_tensor(np.asarray(x), dtype=torch.float32, device=device) - - -def _stack_bags(bags, device): - """ - :param bags: an iterable of (X_bag, prevalence) pairs, all X_bag with the same number of instances - :return: a pair of tensors (X, P) of shape (n_bags, bag_size, n_features) and (n_bags, n_classes) - """ - Xs, ps = zip(*bags) - X = torch.stack([_to_tensor(x, device) for x in Xs]) - P = torch.stack([_to_tensor(p, device) for p in ps]) - return X, P - - -def _mix_two_bags(bag_a, bag_b, bag_size, rng): - """Synthesizes a new bag of size `bag_size` by mixing two given bags with a random ratio, following - the "mixer" idea from the original HistNetQ repo (`UnlabeledMixerBagGenerator`): useful when the - only available training material is a modest number of pre-built samples (e.g., LeQua's dev - samples) and one wants extra intermediate-prevalence bags without access to instance-level labels. - """ - Xa, pa = bag_a - Xb, pb = bag_b - m = rng.random() - na = round(m * bag_size) - nb = bag_size - na - idx_a = rng.choices(range(len(Xa)), k=na) if na > 0 else [] - idx_b = rng.choices(range(len(Xb)), k=nb) if nb > 0 else [] - Xa, Xb = np.asarray(Xa), np.asarray(Xb) - X_mixed = np.concatenate([Xa[idx_a], Xb[idx_b]], axis=0) - p_mixed = m * np.asarray(pa, dtype=float) + (1 - m) * np.asarray(pb, dtype=float) - return X_mixed, p_mixed - - -class HistNetQ(BaseQuantifier): +class HistNetQ(BagTrainedQuantifier): """ Implementation of `HistNetQ `_, a neural network for quantification that learns a differentiable histogram-based representation of a sample, trained @@ -193,14 +113,9 @@ class HistNetQ(BaseQuantifier): end-to-end on whole samples rather than on individually labeled instances, following a symmetric problem setting (learning from bags, predicting on bags). - Training data can be provided in two ways: - - * via :meth:`fit`, from a plain labelled collection (`X`, `y`): training/validation bags are then - generated by resampling from it using a QuaPy sampling protocol (:class:`quapy.protocol.UPP` by - default). - * via :meth:`fit_from_samples`, from a :class:`quapy.protocol.AbstractProtocol` that already yields - the training bags (e.g., :class:`quapy.data._lequa.SamplesFromDir` for LeQua-style pre-built - samples), optionally enriched with synthetic bags mixed from the given ones. + Training data can be provided in two ways: via :meth:`fit`, from a plain labelled collection; or via + :meth:`fit_from_samples`, from a protocol that already yields the training bags. See + :class:`quapy.method._neural_bags.BagTrainedQuantifier` for details. :param feature_extraction_module: a `torch.nn.Module` exposing an `output_size` attribute, used to embed each instance before computing the histogram (e.g., a small MLP for tabular data, a CNN @@ -270,213 +185,37 @@ class HistNetQ(BaseQuantifier): checkpointdir='../checkpoint', checkpointname=None, verbose=0): - self.feature_extraction_module = feature_extraction_module + super().__init__( + feature_extraction_module=feature_extraction_module, + linear_sizes=linear_sizes, + dropout=dropout, + output_function=output_function, + bag_size=bag_size, + n_bags_train=n_bags_train, + n_bags_val=n_bags_val, + train_epochs=train_epochs, + patience=patience, + start_lr=start_lr, + end_lr=end_lr, + lr_factor=lr_factor, + weight_decay=weight_decay, + quant_loss=quant_loss, + batch_size=batch_size, + protocol=protocol, + protocol_params=protocol_params, + val_split=val_split, + device=device, + random_state=random_state, + checkpointdir=checkpointdir, + checkpointname=checkpointname, + verbose=verbose, + ) self.n_bins = n_bins self.quantiles = quantiles - self.linear_sizes = linear_sizes - self.dropout = dropout - self.output_function = output_function - self.bag_size = bag_size - self.n_bags_train = n_bags_train - self.n_bags_val = n_bags_val - self.train_epochs = train_epochs - self.patience = patience - self.start_lr = start_lr - self.end_lr = end_lr - self.lr_factor = lr_factor - self.weight_decay = weight_decay - self.quant_loss = quant_loss if quant_loss is not None else torch.nn.L1Loss() - self.batch_size = batch_size - self.protocol = protocol - self.protocol_params = protocol_params - self.val_split = val_split - self.device = torch.device(device) - self.random_state = random_state - if checkpointname is None: - local_random = random.Random() - random_code = '-'.join(str(local_random.randint(0, 1000000)) for _ in range(5)) - checkpointname = 'HistNetQ-' + random_code - self.checkpointdir = checkpointdir - self.checkpoint = os.path.join(checkpointdir, checkpointname) - self.verbose = verbose - self._classes_ = None @property - def classes_(self): - return self._classes_ + def _checkpoint_prefix(self): + return 'HistNetQ' - def fit(self, X, y): - """ - Trains HistNetQ from a plain labelled collection, generating training and validation bags by - resampling from it via `self.protocol` (a fresh random sequence of bags every epoch for - training, and a fixed, reproducible sequence for validation). - - :param X: the training instances - :param y: the labels of X - :return: self - """ - data = LabelledCollection(X, y) - self._classes_ = data.classes_ - train_data, val_data = data.split_stratified(train_prop=1 - self.val_split, random_state=self.random_state) - - protocol_params = self.protocol_params or {} - - def train_bags(): - sampler = self.protocol( - train_data, sample_size=self.bag_size, repeats=self.n_bags_train, random_state=None, - **protocol_params - ) - return sampler() - - def val_bags(): - sampler = self.protocol( - val_data, sample_size=self.bag_size, repeats=self.n_bags_val, random_state=self.random_state, - **protocol_params - ) - return sampler() - - n_features = train_data.instances.shape[1] - self._fit_loop(train_bags, val_bags, n_features, n_bags_train=self.n_bags_train, n_bags_val=self.n_bags_val) - return self - - def fit_from_samples(self, protocol: AbstractProtocol, val_protocol: AbstractProtocol = None, - mix_bags=False, mix_bags_proportion=0.5): - """ - Trains HistNetQ from a protocol that already yields the training bags (e.g., - :class:`quapy.data._lequa.SamplesFromDir`, for LeQua-style pre-built samples), instead of - resampling from a labelled collection. This is the entry point to use whenever only bags of - known prevalence are available (no instance-level labels). - - :param protocol: an :class:`AbstractProtocol` yielding `(sample, prevalence)` pairs; consumed - once and kept in memory (expected to be of modest size, as is typical of pre-built sample - collections). - :param val_protocol: an optional, separate protocol providing the validation bags; if None, a - `val_split` fraction of the bags returned by `protocol` is held out instead. - :param mix_bags: if True, in addition to the bags returned by `protocol`, synthesize extra bags - each epoch by mixing random pairs of the given bags with a random ratio (a substitute for - the original repo's `UnlabeledMixerBagGenerator`, useful to broaden the coverage of - prevalence values beyond what the given bags exhibit). - :param mix_bags_proportion: proportion (relative to the number of base training bags) of extra - mixed bags to generate per epoch when `mix_bags=True` (default 0.5). - :return: self - """ - assert isinstance(protocol, AbstractProtocol), 'protocol must be an instance of AbstractProtocol' - base_bags = list(protocol()) - n_classes = len(np.asarray(base_bags[0][1])) - self._classes_ = np.arange(n_classes) - - if val_protocol is not None: - val_bags_list = list(val_protocol()) - else: - n_val = max(1, int(len(base_bags) * self.val_split)) - val_bags_list = base_bags[:n_val] - base_bags = base_bags[n_val:] - - rng = random.Random(self.random_state) - n_mixed = round(len(base_bags) * mix_bags_proportion) if mix_bags else 0 - - def train_bags(): - bags = list(base_bags) - if n_mixed > 0: - for _ in range(n_mixed): - a, b = rng.choice(base_bags), rng.choice(base_bags) - bags.append(_mix_two_bags(a, b, self.bag_size, rng)) - rng.shuffle(bags) - return bags - - def val_bags(): - return val_bags_list - - n_features = np.asarray(base_bags[0][0]).shape[1] - self._fit_loop( - train_bags, val_bags, n_features, - n_bags_train=len(base_bags) + n_mixed, n_bags_val=len(val_bags_list) - ) - return self - - def _fit_loop(self, train_bags_fn, val_bags_fn, n_features, n_bags_train, n_bags_val): - os.makedirs(self.checkpointdir, exist_ok=True) - n_classes = len(self._classes_) - - fe = self.feature_extraction_module - if fe is None: - fe = _IdentityFeatureExtractionModule(n_features) - self.histnet = _HistNetModule( - fe, n_classes, n_bins=self.n_bins, quantiles=self.quantiles, linear_sizes=self.linear_sizes, - dropout=self.dropout, output_function=self.output_function - ).to(self.device) - - optim = torch.optim.Adam(self.histnet.parameters(), lr=self.start_lr, weight_decay=self.weight_decay) - scheduler = torch.optim.lr_scheduler.ReduceLROnPlateau(optim, patience=self.patience, factor=self.lr_factor) - early_stop = EarlyStop(self.patience, lower_is_better=True) - - best_state = copy.deepcopy(self.histnet.state_dict()) - for epoch in range(self.train_epochs): - tr_loss = self._run_epoch(train_bags_fn(), n_bags_train, optim, train=True, epoch=epoch) - va_loss = self._run_epoch(val_bags_fn(), n_bags_val, optim=None, train=False, epoch=epoch) - - early_stop(va_loss, epoch) - if early_stop.IMPROVED: - best_state = copy.deepcopy(self.histnet.state_dict()) - torch.save(best_state, self.checkpoint) - elif early_stop.STOP: - if self.verbose > 0: - print(f'[HistNetQ] training ended by patience exhausted at epoch {epoch}; ' - f'restoring best model from epoch {early_stop.best_epoch}') - break - - scheduler.step(va_loss) - if optim.param_groups[0]['lr'] < self.end_lr: - if self.verbose > 0: - print(f'[HistNetQ] early stopping in epoch {epoch} (learning rate below end_lr)') - break - - self.histnet.load_state_dict(best_state) - - def _run_epoch(self, bags, n_bags, optim, train, epoch): - self.histnet.train(mode=train) - losses = [] - pbar = tqdm(bags, total=n_bags, disable=self.verbose == 0) - batch = [] - - def process_batch(batch): - X, P = _stack_bags(batch, self.device) - if train: - optim.zero_grad() - P_hat = self.histnet.forward(X) - loss = self.quant_loss(P, P_hat) - loss.backward() - optim.step() - else: - with torch.no_grad(): - P_hat = self.histnet.forward(X) - loss = self.quant_loss(P, P_hat) - return loss.item() - - for bag in pbar: - batch.append(bag) - if len(batch) == self.batch_size: - losses.append(process_batch(batch)) - batch = [] - pbar.set_description( - f'[HistNetQ] epoch={epoch} {"train" if train else "val"}-loss={np.mean(losses):.5f}' - ) - if batch: - losses.append(process_batch(batch)) - - return np.mean(losses) if losses else float('inf') - - def predict(self, X): - """ - Generates a class prevalence estimate for the sample `X`, via a single forward pass of the - trained network (the histogram layer aggregates over however many instances are given, so `X` - need not match the `bag_size` used during training). - - :param X: the test instances - :return: `np.ndarray` of shape `(n_classes,)` with the class prevalence estimates - """ - self.histnet.eval() - with torch.no_grad(): - X_t = _to_tensor(X, self.device).unsqueeze(0) - prevalence = self.histnet.forward(X_t) - return prevalence.cpu().numpy().flatten() + def _build_quantmodule(self, n_features): + return _SigmoidHistogram(n_features=n_features, n_bins=self.n_bins, quantiles=self.quantiles) diff --git a/quapy/method/_neural_bags.py b/quapy/method/_neural_bags.py new file mode 100644 index 0000000..f32af71 --- /dev/null +++ b/quapy/method/_neural_bags.py @@ -0,0 +1,443 @@ +""" +Shared machinery for QuaPy's "bag-trained" neural quantifiers, i.e., methods that -- like +:class:`quapy.method.meta.QuaNet` -- do not follow the classify-then-aggregate pattern of +:class:`quapy.method.aggregative.AggregativeQuantifier`, but are instead trained and evaluated +end-to-end on whole samples ("bags") of known prevalence. + +:class:`quapy.method._histnet.HistNetQ` and :class:`quapy.method._gmnet.GMNet` share the same overall +architecture (feature_extraction -> quantification module -> small MLP head -> softmax/normalize) and +the same bag-based training protocol (bag generation via a QuaPy sampling protocol, early stopping, LR +scheduling, checkpointing). This module factors that common part out into :class:`BagTrainedQuantifier`; +concrete subclasses only need to supply the quantification module placed between the feature extractor +and the shared head (see :meth:`BagTrainedQuantifier._build_quantmodule`). +""" +import copy +import os +import random +from abc import ABC, abstractmethod + +import numpy as np +import torch +import torch.nn as nn +from tqdm import tqdm + +from quapy.data import LabelledCollection +from quapy.method.base import BaseQuantifier +from quapy.protocol import AbstractProtocol, UPP +from quapy.util import EarlyStop + + +class IdentityFeatureExtractionModule(nn.Module): + """Used when no feature extraction module is provided: instances are assumed to already be in + their final numeric representation.""" + + def __init__(self, input_size): + super().__init__() + self.output_size = input_size + + def forward(self, x): + return x + + +def to_tensor(x, device): + if torch.is_tensor(x): + return x.to(device=device, dtype=torch.float32) + if hasattr(x, 'toarray'): # scipy sparse + x = x.toarray() + return torch.as_tensor(np.asarray(x), dtype=torch.float32, device=device) + + +def stack_bags(bags, device): + """ + :param bags: an iterable of (X_bag, prevalence) pairs, all X_bag with the same number of instances + :return: a pair of tensors (X, P) of shape (n_bags, bag_size, n_features) and (n_bags, n_classes) + """ + Xs, ps = zip(*bags) + X = torch.stack([to_tensor(x, device) for x in Xs]) + P = torch.stack([to_tensor(p, device) for p in ps]) + return X, P + + +def mix_two_bags(bag_a, bag_b, bag_size, rng): + """Synthesizes a new bag of size `bag_size` by mixing two given bags with a random ratio, following + the "mixer" idea from the original HistNetQ repo (`UnlabeledMixerBagGenerator`): useful when the + only available training material is a modest number of pre-built samples (e.g., LeQua's dev + samples) and one wants extra intermediate-prevalence bags without access to instance-level labels. + """ + Xa, pa = bag_a + Xb, pb = bag_b + m = rng.random() + na = round(m * bag_size) + nb = bag_size - na + idx_a = rng.choices(range(len(Xa)), k=na) if na > 0 else [] + idx_b = rng.choices(range(len(Xb)), k=nb) if nb > 0 else [] + Xa, Xb = np.asarray(Xa), np.asarray(Xb) + X_mixed = np.concatenate([Xa[idx_a], Xb[idx_b]], axis=0) + p_mixed = m * np.asarray(pa, dtype=float) + (1 - m) * np.asarray(pb, dtype=float) + return X_mixed, p_mixed + + +def build_output_head(input_size, n_classes, linear_sizes, dropout, output_function): + """Builds the small MLP + output activation shared by every bag-trained quantifier's head: a stack + of (Linear, LeakyReLU, Dropout) blocks sized by `linear_sizes`, followed by a final Linear to + `n_classes` and either a softmax or an L1-normalization (applied in :class:`BagNetworkModule`), both + yielding a valid prevalence vector.""" + output_module = nn.Sequential() + prev_size = input_size + for i, linear_size in enumerate(linear_sizes): + output_module.add_module(f'linear_{i}', nn.Linear(prev_size, linear_size)) + output_module.add_module(f'leakyrelu_{i}', nn.LeakyReLU()) + output_module.add_module(f'dropout_{i}', nn.Dropout(dropout)) + prev_size = linear_size + output_module.add_module('last_linear', nn.Linear(prev_size, n_classes)) + if output_function == 'softmax': + output_module.add_module('softmax', nn.Softmax(dim=1)) + elif output_function == 'normalize': + output_module.add_module('relu', nn.ReLU()) + else: + raise ValueError(f"unknown {output_function=}; valid ones are 'softmax', 'normalize'") + return output_module + + +class BagNetworkModule(nn.Module): + """The full network shared by every bag-trained quantifier: feature extraction, a pluggable + quantification module (mapping a bag of instance-level features to a single per-bag + representation), and the shared MLP head. + + :param quantmodule: a `torch.nn.Module` exposing an `output_size` attribute, mapping a tensor of + shape (batch_size, bag_size, n_features) to one of shape (batch_size, quantmodule.output_size). + """ + + def __init__(self, feature_extraction_module, quantmodule, n_classes, linear_sizes, dropout, output_function): + super().__init__() + self.feature_extraction_module = feature_extraction_module + self.quantmodule = quantmodule + self.output_function = output_function + self.output_module = build_output_head( + quantmodule.output_size, n_classes, linear_sizes, dropout, output_function + ) + + def forward(self, bag): + # bag: (batch_size, bag_size, n_features) + features = self.feature_extraction_module(bag) + representation = self.quantmodule(features) + out = self.output_module(representation) + if self.output_function == 'normalize': + out = nn.functional.normalize(out, p=1, dim=1) + return out + + +class BagTrainedQuantifier(BaseQuantifier, ABC): + """ + Base class for QuaPy's neural quantifiers trained end-to-end on whole samples ("bags") of known + prevalence, rather than following the classify-then-aggregate pattern of + :class:`quapy.method.aggregative.AggregativeQuantifier` (in the spirit of + :class:`quapy.method.meta.QuaNet`). Concrete subclasses only need to provide the quantification + module placed between the feature extractor and the shared MLP head (see + :meth:`_build_quantmodule`) and a checkpoint-name prefix (see :attr:`_checkpoint_prefix`); bag + generation, the training/validation loop, early stopping, LR scheduling, checkpointing, and + prediction are all shared. + + Training data can be provided in two ways: + + * via :meth:`fit`, from a plain labelled collection (`X`, `y`): training/validation bags are then + generated by resampling from it using a QuaPy sampling protocol (:class:`quapy.protocol.UPP` by + default). + * via :meth:`fit_from_samples`, from a :class:`quapy.protocol.AbstractProtocol` that already yields + the training bags (e.g., :class:`quapy.data._lequa.SamplesFromDir` for LeQua-style pre-built + samples), optionally enriched with synthetic bags mixed from the given ones. + + :param feature_extraction_module: a `torch.nn.Module` exposing an `output_size` attribute, used to + embed each instance before the quantification module (e.g., a small MLP for tabular data, a CNN + for images). If None (default), an identity module is used, i.e., the instances in `X` are + assumed to already be in their final numeric representation. + :param linear_sizes: tuple of ints with the sizes of the linear layers used in the shared head + (default empty, i.e., only the final classification layer is used). + :param dropout: dropout applied after each of the `linear_sizes` layers (default 0). + :param output_function: either 'softmax' or 'normalize' (L1); both yield a valid prevalence vector + (default 'softmax'). + :param bag_size: number of instances per training/validation bag (default 500). + :param n_bags_train: number of bags generated per training epoch (default 500). + :param n_bags_val: number of bags generated per validation epoch (default 500). + :param train_epochs: maximum number of training epochs (default 200). + :param patience: number of epochs without improvement in validation loss before early-stopping + (default 20). + :param start_lr: initial learning rate (default 1e-3). + :param end_lr: once the learning rate decays below this value, training stops (default 1e-6). + :param lr_factor: factor by which the learning rate is reduced after `patience` epochs without + improvement (default 0.1). + :param weight_decay: L2 regularization (default 0). + :param quant_loss: the quantification loss to minimize (default `torch.nn.L1Loss()`), called as + `quant_loss(true_prevalences, predicted_prevalences)`. + :param batch_size: number of bags per gradient update (default 16). + :param protocol: the :class:`quapy.protocol.AbstractStochasticSeededProtocol` subclass used by + :meth:`fit` to resample bags from the given labelled collection (default + :class:`quapy.protocol.UPP`, which draws bags with prevalence sampled uniformly at random from + the simplex). + :param protocol_params: dict of extra keyword arguments passed to `protocol` (besides `data`, + `sample_size`, `repeats`, and `random_state`, which are set internally); default None. + :param val_split: float in (0,1), the proportion of the collection given to :meth:`fit` that is held + out (via stratified sampling) for validation and early stopping (default 0.4). + :param device: `'cpu'` or `'cuda'` (default 'cpu'). + :param random_state: seed used for the train/validation split and for the (fixed) validation + sampling sequence (default 0). + :param checkpointdir: directory where the best model found during training is stored (default + '../checkpoint'). + :param checkpointname: name of the checkpoint file; if None (default), a random name prefixed by + :attr:`_checkpoint_prefix` is generated. + :param verbose: verbosity level; if >0, shows a progress bar with the current losses (default 0). + """ + + def __init__(self, + feature_extraction_module=None, + linear_sizes=(), + dropout=0., + output_function='softmax', + bag_size=500, + n_bags_train=500, + n_bags_val=500, + train_epochs=200, + patience=20, + start_lr=1e-3, + end_lr=1e-6, + lr_factor=0.1, + weight_decay=0., + quant_loss=None, + batch_size=16, + protocol=UPP, + protocol_params=None, + val_split=0.4, + device='cpu', + random_state=0, + checkpointdir='../checkpoint', + checkpointname=None, + verbose=0): + self.feature_extraction_module = feature_extraction_module + self.linear_sizes = linear_sizes + self.dropout = dropout + self.output_function = output_function + self.bag_size = bag_size + self.n_bags_train = n_bags_train + self.n_bags_val = n_bags_val + self.train_epochs = train_epochs + self.patience = patience + self.start_lr = start_lr + self.end_lr = end_lr + self.lr_factor = lr_factor + self.weight_decay = weight_decay + self.quant_loss = quant_loss if quant_loss is not None else torch.nn.L1Loss() + self.batch_size = batch_size + self.protocol = protocol + self.protocol_params = protocol_params + self.val_split = val_split + self.device = torch.device(device) + self.random_state = random_state + if checkpointname is None: + local_random = random.Random() + random_code = '-'.join(str(local_random.randint(0, 1000000)) for _ in range(5)) + checkpointname = f'{self._checkpoint_prefix}-{random_code}' + self.checkpointdir = checkpointdir + self.checkpoint = os.path.join(checkpointdir, checkpointname) + self.verbose = verbose + self._classes_ = None + self.model = None + + @property + def classes_(self): + return self._classes_ + + @property + @abstractmethod + def _checkpoint_prefix(self): + """Short name used as the default checkpoint filename prefix (e.g. 'HistNetQ', 'GMNet').""" + ... + + @abstractmethod + def _build_quantmodule(self, n_features): + """Builds the module placed between the (already feature-extracted) instances and the shared + MLP head. Must expose an `output_size` attribute and accept input of shape + (batch_size, bag_size, n_features), returning one of shape (batch_size, output_size).""" + ... + + def _extra_loss(self): + """Optional additional term added to the quantification loss during training (e.g., GMNet's CKA + regularization across GM layers). Returns 0 by default.""" + return 0. + + def fit(self, X, y): + """ + Trains the quantifier from a plain labelled collection, generating training and validation bags + by resampling from it via `self.protocol` (a fresh random sequence of bags every epoch for + training, and a fixed, reproducible sequence for validation). + + :param X: the training instances + :param y: the labels of X + :return: self + """ + data = LabelledCollection(X, y) + self._classes_ = data.classes_ + train_data, val_data = data.split_stratified(train_prop=1 - self.val_split, random_state=self.random_state) + + protocol_params = self.protocol_params or {} + + def train_bags(): + sampler = self.protocol( + train_data, sample_size=self.bag_size, repeats=self.n_bags_train, random_state=None, + **protocol_params + ) + return sampler() + + def val_bags(): + sampler = self.protocol( + val_data, sample_size=self.bag_size, repeats=self.n_bags_val, random_state=self.random_state, + **protocol_params + ) + return sampler() + + n_features = train_data.instances.shape[1] + self._fit_loop(train_bags, val_bags, n_features, n_bags_train=self.n_bags_train, n_bags_val=self.n_bags_val) + return self + + def fit_from_samples(self, protocol: AbstractProtocol, val_protocol: AbstractProtocol = None, + mix_bags=False, mix_bags_proportion=0.5): + """ + Trains the quantifier from a protocol that already yields the training bags (e.g., + :class:`quapy.data._lequa.SamplesFromDir`, for LeQua-style pre-built samples), instead of + resampling from a labelled collection. This is the entry point to use whenever only bags of + known prevalence are available (no instance-level labels). + + :param protocol: an :class:`AbstractProtocol` yielding `(sample, prevalence)` pairs; consumed + once and kept in memory (expected to be of modest size, as is typical of pre-built sample + collections). + :param val_protocol: an optional, separate protocol providing the validation bags; if None, a + `val_split` fraction of the bags returned by `protocol` is held out instead. + :param mix_bags: if True, in addition to the bags returned by `protocol`, synthesize extra bags + each epoch by mixing random pairs of the given bags with a random ratio (a substitute for + the original HistNetQ repo's `UnlabeledMixerBagGenerator`, useful to broaden the coverage of + prevalence values beyond what the given bags exhibit). + :param mix_bags_proportion: proportion (relative to the number of base training bags) of extra + mixed bags to generate per epoch when `mix_bags=True` (default 0.5). + :return: self + """ + assert isinstance(protocol, AbstractProtocol), 'protocol must be an instance of AbstractProtocol' + base_bags = list(protocol()) + n_classes = len(np.asarray(base_bags[0][1])) + self._classes_ = np.arange(n_classes) + + if val_protocol is not None: + val_bags_list = list(val_protocol()) + else: + n_val = max(1, int(len(base_bags) * self.val_split)) + val_bags_list = base_bags[:n_val] + base_bags = base_bags[n_val:] + + rng = random.Random(self.random_state) + n_mixed = round(len(base_bags) * mix_bags_proportion) if mix_bags else 0 + + def train_bags(): + bags = list(base_bags) + if n_mixed > 0: + for _ in range(n_mixed): + a, b = rng.choice(base_bags), rng.choice(base_bags) + bags.append(mix_two_bags(a, b, self.bag_size, rng)) + rng.shuffle(bags) + return bags + + def val_bags(): + return val_bags_list + + n_features = np.asarray(base_bags[0][0]).shape[1] + self._fit_loop( + train_bags, val_bags, n_features, + n_bags_train=len(base_bags) + n_mixed, n_bags_val=len(val_bags_list) + ) + return self + + def _fit_loop(self, train_bags_fn, val_bags_fn, n_features, n_bags_train, n_bags_val): + os.makedirs(self.checkpointdir, exist_ok=True) + n_classes = len(self._classes_) + + fe = self.feature_extraction_module + if fe is None: + fe = IdentityFeatureExtractionModule(n_features) + quantmodule = self._build_quantmodule(fe.output_size) + self.model = BagNetworkModule( + fe, quantmodule, n_classes, linear_sizes=self.linear_sizes, dropout=self.dropout, + output_function=self.output_function + ).to(self.device) + + optim = torch.optim.Adam(self.model.parameters(), lr=self.start_lr, weight_decay=self.weight_decay) + scheduler = torch.optim.lr_scheduler.ReduceLROnPlateau(optim, patience=self.patience, factor=self.lr_factor) + early_stop = EarlyStop(self.patience, lower_is_better=True) + + best_state = copy.deepcopy(self.model.state_dict()) + for epoch in range(self.train_epochs): + self._run_epoch(train_bags_fn(), n_bags_train, optim, train=True, epoch=epoch) + va_loss = self._run_epoch(val_bags_fn(), n_bags_val, optim=None, train=False, epoch=epoch) + + early_stop(va_loss, epoch) + if early_stop.IMPROVED: + best_state = copy.deepcopy(self.model.state_dict()) + torch.save(best_state, self.checkpoint) + elif early_stop.STOP: + if self.verbose > 0: + print(f'[{self._checkpoint_prefix}] training ended by patience exhausted at epoch {epoch}; ' + f'restoring best model from epoch {early_stop.best_epoch}') + break + + scheduler.step(va_loss) + if optim.param_groups[0]['lr'] < self.end_lr: + if self.verbose > 0: + print(f'[{self._checkpoint_prefix}] early stopping in epoch {epoch} (learning rate below end_lr)') + break + + self.model.load_state_dict(best_state) + + def _run_epoch(self, bags, n_bags, optim, train, epoch): + self.model.train(mode=train) + losses = [] + pbar = tqdm(bags, total=n_bags, disable=self.verbose == 0) + batch = [] + + def process_batch(batch): + X, P = stack_bags(batch, self.device) + if train: + optim.zero_grad() + P_hat = self.model.forward(X) + loss = self.quant_loss(P, P_hat) + self._extra_loss() + loss.backward() + optim.step() + else: + with torch.no_grad(): + P_hat = self.model.forward(X) + loss = self.quant_loss(P, P_hat) + return loss.item() + + for bag in pbar: + batch.append(bag) + if len(batch) == self.batch_size: + losses.append(process_batch(batch)) + batch = [] + pbar.set_description( + f'[{self._checkpoint_prefix}] epoch={epoch} {"train" if train else "val"}-' + f'loss={np.mean(losses):.5f}' + ) + if batch: + losses.append(process_batch(batch)) + + return np.mean(losses) if losses else float('inf') + + def predict(self, X): + """ + Generates a class prevalence estimate for the sample `X`, via a single forward pass of the + trained network (the quantification module aggregates over however many instances are given, so + `X` need not match the `bag_size` used during training). + + :param X: the test instances + :return: `np.ndarray` of shape `(n_classes,)` with the class prevalence estimates + """ + self.model.eval() + with torch.no_grad(): + X_t = to_tensor(X, self.device).unsqueeze(0) + prevalence = self.model.forward(X_t) + return prevalence.cpu().numpy().flatten() diff --git a/quapy/method/meta.py b/quapy/method/meta.py index 7d15b7c..58c26dd 100644 --- a/quapy/method/meta.py +++ b/quapy/method/meta.py @@ -37,6 +37,17 @@ if _histnet: else: HistNetQ = "HistNetQ is not available due to missing torch package" +try: + from . import _gmnet +except ModuleNotFoundError: + _gmnet = None + + +if _gmnet: + GMNet = _gmnet.GMNet +else: + GMNet = "GMNet is not available due to missing torch and/or geotorch packages" + class MedianEstimator(BinaryQuantifier): """ diff --git a/quapy/tests/test_methods.py b/quapy/tests/test_methods.py index 7a225dd..0e41ff5 100644 --- a/quapy/tests/test_methods.py +++ b/quapy/tests/test_methods.py @@ -127,6 +127,38 @@ class TestMethods(unittest.TestCase): estim_prevalences2 = model2.predict(dataset.test.X) self.assertTrue(check_prevalence_vector(estim_prevalences2)) + def test_gmnet(self): + try: + import torch + import geotorch + except ModuleNotFoundError: + print('the torch and/or geotorch packages are not installed; skipping unit test for GMNet') + return + + from quapy.method.meta import GMNet + from quapy.protocol import UPP + + for dataset in TestMethods.datasets: + # single GM layer, no CKA regularization + model = GMNet( + bag_size=20, n_bags_train=10, n_bags_val=5, train_epochs=2, patience=1, batch_size=2, + device='cpu', checkpointdir='./checkpoint_test_gmnet' + ) + model.fit(*dataset.training.Xy) + estim_prevalences = model.predict(dataset.test.X) + self.assertTrue(check_prevalence_vector(estim_prevalences)) + + # multiple GM layers + CKA regularization, and fit_from_samples + given_samples = UPP(dataset.training, sample_size=20, repeats=8, random_state=1) + val_samples = UPP(dataset.training, sample_size=20, repeats=4, random_state=2) + model2 = GMNet( + n_gm_layers=2, num_gaussians=3, gaussian_dimensions=4, cka_regularization=0.1, + bag_size=20, train_epochs=2, patience=1, batch_size=2, device='cpu', + checkpointdir='./checkpoint_test_gmnet' + ) + model2.fit_from_samples(given_samples, val_protocol=val_samples, mix_bags=True) + estim_prevalences2 = model2.predict(dataset.test.X) + self.assertTrue(check_prevalence_vector(estim_prevalences2)) def test_composable(self): try: diff --git a/setup.py b/setup.py index 22c424f..3ff18da 100644 --- a/setup.py +++ b/setup.py @@ -131,7 +131,7 @@ setup( # projects. extras_require={ # Optional 'bayes': ['jax', 'jaxlib', 'numpyro', 'pystan', 'setuptools<82'], - 'neural': ['torch'], + 'neural': ['torch', 'geotorch'], 'tests': ['certifi'], 'docs' : ['pydata-sphinx-theme', 'myst-parser', 'sphinx-design'], }, From 2a50fe902c7aaf47f0a2aad991ab51c6e7a079a7 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Tue, 8 Sep 2026 14:47:13 +0200 Subject: [PATCH 4/9] Fix broken imports left over from the BBSE/LEIP commit quapy/classification/__init__.py imported a labelshift module that doesn't exist anywhere in the repo (and nothing else references it), breaking `import quapy` entirely; the import is removed. quapy/method/aggregative.py imported a nonexistent `_liep` module instead of the actual file, _liep_draft.py, which itself never got a LEIP class added (only helper functions such as leip()). Guards the import so LEIP degrades to an "not available" placeholder, like the other optional neural methods, instead of crashing the whole package. test_leip/test_leip_fixed_tau still fail as a result; finishing LEIP is left for a follow-up. Co-Authored-By: Claude Sonnet 5 --- quapy/classification/__init__.py | 1 - quapy/method/aggregative.py | 8 +++++--- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/quapy/classification/__init__.py b/quapy/classification/__init__.py index 33551b7..b63e0eb 100644 --- a/quapy/classification/__init__.py +++ b/quapy/classification/__init__.py @@ -1,4 +1,3 @@ from . import calibration -from . import labelshift from . import methods from . import svmperf diff --git a/quapy/method/aggregative.py b/quapy/method/aggregative.py index 23341b6..1decb39 100644 --- a/quapy/method/aggregative.py +++ b/quapy/method/aggregative.py @@ -2218,9 +2218,11 @@ KDEyML = _kdey.KDEyML KDEyHD = _kdey.KDEyHD KDEyCS = _kdey.KDEyCS -from . import _liep - -LEIP = _liep.LEIP +try: + from . import _liep_draft as _liep + LEIP = _liep.LEIP +except AttributeError: + LEIP = "LEIP is not available (incomplete implementation in _liep_draft.py)" # --------------------------------------------------------------- From cb73c78cf3509c6db5c435cc2524518cd10c8f5e Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Wed, 9 Sep 2026 11:40:12 +0200 Subject: [PATCH 5/9] Add manual entry for GMNet Documents GMNet in docs/source/manuals/methods.md, mirroring the existing HistNetQ entry: paper reference, torch/geotorch requirement, a basic usage example, and a second example showing multi-branch + CKA regularization. Co-Authored-By: Claude Sonnet 5 --- docs/source/manuals/methods.md | 40 ++++++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) diff --git a/docs/source/manuals/methods.md b/docs/source/manuals/methods.md index 5be168f..26fea42 100644 --- a/docs/source/manuals/methods.md +++ b/docs/source/manuals/methods.md @@ -934,6 +934,46 @@ HistNetQ can alternatively be trained directly from a protocol that already prov samples (e.g., when only bag-level prevalence values are available), via the `fit_from_samples` method; see the API documentation for further details. +### GMNet + +QuaPy offers an implementation of GMNet, a deep learning model that represents each instance of a +bag by its likelihood under one or more learned mixtures of Gaussians, presented in: + +[_Pérez-Mon, O., del Coz, J.J., & González, P. (2026). +Quantification via Gaussian latent space representations. +Neural Networks._](https://www.sciencedirect.com/science/article/pii/S0893608026003473) + +This model requires `torch` and `geotorch` to be installed (`geotorch` is used to keep the Gaussian +layers' covariance matrices positive-definite during training). Like HistNetQ, GMNet is trained +end-to-end on samples ("bags") of known prevalence rather than on individually labeled instances, +and requires no classifier, only an optional feature extraction module. + +```python +import quapy as qp +from quapy.method.meta import GMNet + +dataset = qp.datasets.fetch_UCIBinaryDataset('haberman') + +model = GMNet(bag_size=100, device='cpu') +model.fit(*dataset.training.Xy) +estim_prevalence = model.predict(dataset.test.X) +``` + +GMNet supports stacking multiple "GM branches" (via `n_gm_layers`, `num_gaussians`, and +`gaussian_dimensions`), optionally regularized with a CKA (Centered Kernel Alignment) term that +encourages the branches to learn dissimilar latent representations, set through +`cka_regularization`: + +```python +model = GMNet( + n_gm_layers=2, num_gaussians=(4, 4), gaussian_dimensions=(8, 8), + cka_regularization=0.1, bag_size=100, device='cpu' +) +``` + +Like HistNetQ, GMNet can alternatively be trained directly from a protocol that already provides the +training samples, via the `fit_from_samples` method; see the API documentation for further details. + ## Quantifiers with Uncertainty Quantification From 3ecbadb921419d3287fefc3bc95354ef9090dda1 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Fri, 18 Sep 2026 11:54:50 +0200 Subject: [PATCH 6/9] adding liep --- docs/source/manuals/label-shift-adaptation.md | 108 +++++++++ quapy/classification/labelshift.py | 110 +++++++++ quapy/method/_liep.py | 227 ++++++++++++++++++ quapy/tests/test_labelshift.py | 84 +++++++ 4 files changed, 529 insertions(+) create mode 100644 docs/source/manuals/label-shift-adaptation.md create mode 100644 quapy/classification/labelshift.py create mode 100644 quapy/method/_liep.py create mode 100644 quapy/tests/test_labelshift.py diff --git a/docs/source/manuals/label-shift-adaptation.md b/docs/source/manuals/label-shift-adaptation.md new file mode 100644 index 0000000..f71994f --- /dev/null +++ b/docs/source/manuals/label-shift-adaptation.md @@ -0,0 +1,108 @@ +# Label Shift Adaptation + +Quantification methods estimate the class prevalence values of a test sample, but a prevalence estimate is +not, by itself, a classifier that has been corrected for the shift it describes. In some scenarios you may +want the latter: a classifier whose decisions account for the (estimated) change in class priors between +the training and the test distributions, rather than a single prevalence vector. + +QuaPy provides two building blocks for this: + +* `qp.method.aggregative.ImportanceWeightQuantifier`, an interface implemented by the quantifiers that + natively compute a vector of importance weights as part of estimating the target prevalence. +* `qp.classification.labelshift.LabelShiftedClassifier`, a classifier wrapper that uses any quantifier + (whether or not it implements the interface above) to adapt a classifier to a specific, (potentially) + shifted batch of test instances. + +## The `ImportanceWeightQuantifier` interface + +Some quantifiers (currently, `RLLS`, `BBSEhard`, and `BBSEsoft`; see the +{ref}`Regularized Learning under Label Shift (RLLS) ` +and {ref}`Black Box Shift Estimation (BBSE) ` sections of +the methods manual) work by first estimating a vector of importance weights + +:::{math} +w_y = \frac{Q(y)}{P(y)} +::: + +with :math:`P` and :math:`Q` the training and target class distributions, and only then rescaling the +training prevalence by these weights to obtain the target prevalence estimate: +:math:`\hat{Q}(y) \propto w_y \cdot P(y)`. For these methods, the weight vector is not an afterthought +computed from the prevalence estimate; it is the primitive quantity from which the estimate itself is +derived. + +`ImportanceWeightQuantifier` exposes this quantity directly, in addition to the usual `predict` method +inherited from any quantifier: + +* `get_importance_weights(instances)` returns the weight vector $w$ estimated for a given batch of + (unlabelled) target instances. +* `quantify_and_weigh(instances)` returns a tuple `(prevalence, weights)`, computed from a single pass of + classifier predictions over the instances (avoiding the double classification that calling `predict` and + `get_importance_weights` separately would incur). + +Neither method mutates any internal state of the quantifier, so both are safe to call concurrently, on the +same fitted instance, for different batches of target instances. + +```python +import quapy as qp +from quapy.method.aggregative import BBSEhard +from sklearn.linear_model import LogisticRegression + +train, test = qp.datasets.fetch_UCIBinaryDataset('haberman').train_test + +model = BBSEhard(LogisticRegression(max_iter=2000), val_split=5) +model.fit(*train.Xy) + +prevalence, weights = model.quantify_and_weigh(test.X) +``` + +## Adapting a classifier: `LabelShiftedClassifier` + +Rescaling the classifier's posterior probabilities after the fact (as, e.g., `EMQ` does) is one way of +accounting for an estimated prevalence shift. Another common practice is instead to retrain the classifier +with a reweighted loss, via scikit-learn's `class_weight` constructor parameter; this changes the decision +function of the classifier itself, rather than only rescaling its outputs, which can matter for classifiers +whose decision boundary is not a simple function of the posterior probabilities alone. + +`LabelShiftedClassifier` automates this second strategy. Given a base classifier `h` and any `quapy` +quantifier `q`, it fits `q` on the labelled training data (keeping a copy of it), and, for a given batch of +(unlabelled) target instances, estimates the importance weights: directly, via `get_importance_weights`, if +`q` is an `ImportanceWeightQuantifier`; or otherwise by dividing `q`'s estimated target prevalence by the +training prevalence. It then retrains a *fresh copy* of `h`, with `class_weight` set according to these +weights, on the original training data, and hands off classification of the target instances to this +adapted classifier. Neither the classifier instance `h` passed at construction time, nor `q`'s own internal +classifier (if any), are ever mutated. + +```python +import quapy as qp +from quapy.classification.labelshift import LabelShiftedClassifier +from quapy.method.aggregative import ACC +from sklearn.linear_model import LogisticRegression + +train, test = qp.datasets.fetch_UCIBinaryDataset('haberman').train_test + +h = LogisticRegression(max_iter=2000) +q = ACC(LogisticRegression(max_iter=2000), val_split=5) + +adapter = LabelShiftedClassifier(h, q) +adapter.fit(*train.Xy) + +# label predictions and posterior probabilities from a classifier adapted to test's estimated distribution +y_pred = adapter.predict(test.X) +y_proba = adapter.predict_proba(test.X) + +# or, to obtain the adapted classifier itself +adapted_h = adapter.get_classifier(test.X) +``` + +`h` is required to accept a `class_weight` parameter in its constructor, as `LogisticRegression` and many +other scikit-learn classifiers do; `q` can be any (fitted or unfitted) `quapy` quantifier, and need not be an +`ImportanceWeightQuantifier`. + +Since the adaptation is specific to the batch of target instances given, `predict` and `predict_proba` +retrain `h` on every call. If you plan to issue several predictions against the same (suspected) shifted +distribution, call `get_classifier` once and reuse the classifier instance it returns, rather than calling +`predict`/`predict_proba` repeatedly on the same or related samples. + +A small positive floor, controlled by the `weight_epsilon` constructor parameter (default `1e-4`), is +applied to the estimated weights before they are used as `class_weight`; this keeps a class with training +support from being assigned a weight of exactly `0`, which would silently drop it from the retraining loss. diff --git a/quapy/classification/labelshift.py b/quapy/classification/labelshift.py new file mode 100644 index 0000000..6225f1e --- /dev/null +++ b/quapy/classification/labelshift.py @@ -0,0 +1,110 @@ +import numpy as np +from sklearn.base import BaseEstimator, clone + +import quapy.functional as F + + +class LabelShiftedClassifier(BaseEstimator): + """ + Adapts a classifier to a shifted (unlabelled) target distribution by means of a `quapy` quantifier. + + Given a base classifier `h` and a quantifier `q`, this wrapper: (i) at `fit` time, fits `q` on the labelled + training data and keeps a copy of this data; (ii) at prediction time, for a given batch of (unlabelled) + target instances, estimates the importance weights :math:`w_y=Q(y)/P(y)` that would explain the shift + between the training distribution and the target one (directly, if `q` is an + :class:`quapy.method.aggregative.ImportanceWeightQuantifier`, or otherwise by dividing `q`'s prevalence + estimate by the training prevalence); and (iii) retrains a fresh copy of `h`, reweighted according to `w` + via `class_weight`, on the original training data, before delegating the classification of the target + instances to this adapted classifier. + + Note that a fresh copy of `h` is retrained for every batch of target instances passed to `predict`, + `predict_proba`, or `get_classifier`, since the adaptation is specific to that batch's estimated shift; `q` + and its own internal classifier (if any) are otherwise left untouched. `h` itself is required to accept a + `class_weight` parameter in its constructor (e.g., `sklearn.linear_model.LogisticRegression`). + + :param classifier: a scikit-learn classifier that accepts a `class_weight` parameter in its constructor + :param quantifier: a `quapy` quantifier, used to estimate the target prevalence (or, if it is an + :class:`quapy.method.aggregative.ImportanceWeightQuantifier`, the importance weights directly) + :param weight_epsilon: float, a small positive floor applied to the estimated importance weights (default + 1e-4), so that a class with training support is never assigned a `class_weight` of exactly 0 (which + would silently drop it from the retraining loss) + """ + + def __init__(self, classifier, quantifier, weight_epsilon=1e-4): + self.classifier = classifier + self.quantifier = quantifier + self.weight_epsilon = weight_epsilon + + def fit(self, X, y): + """ + Fits the internal quantifier on the training data, and stores the training data for the subsequent, + per-target-batch, classifier adaptation. + + :param X: array-like of shape `(n_samples, n_features)` with the training instances + :param y: array-like of shape `(n_samples,)` with the class labels + :return: self + """ + self.quantifier.fit(X, y) + self.classes_ = self.quantifier.classes_ + self.train_prevalence_ = F.prevalence_from_labels(y, classes=self.classes_) + self.X_, self.y_ = X, y + return self + + def get_importance_weights(self, X): + """ + Estimates the vector of importance weights :math:`w_y=Q(y)/P(y)` that would explain the shift between + the training distribution and the distribution of the given target instances. This method does not + mutate any internal state, so it is safe to call concurrently for different batches of target + instances. + + :param X: array-like of shape `(n_samples, n_features)` with the (unlabelled) target instances + :return: np.ndarray of shape `(n_classes,)` + """ + from quapy.method.aggregative import ImportanceWeightQuantifier + if isinstance(self.quantifier, ImportanceWeightQuantifier): + weights = self.quantifier.get_importance_weights(X) + else: + test_prevalence = self.quantifier.predict(X) + weights = test_prevalence / np.maximum(self.train_prevalence_, self.weight_epsilon) + return np.clip(weights, self.weight_epsilon, None) + + def get_classifier(self, X): + """ + Returns a fresh classifier instance, sharing `classifier`'s hyperparameters, retrained on the original + training data with `class_weight` set according to the importance weights estimated for the given + target instances. The classifier passed at construction time, and the one used internally by the + quantifier, are left untouched. + + :param X: array-like of shape `(n_samples, n_features)` with the (unlabelled) target instances + :return: a fitted scikit-learn classifier + """ + weights = self.get_importance_weights(X) + class_weight = dict(zip(self.classes_, weights)) + adapted = clone(self.classifier) + adapted.set_params(class_weight=class_weight) + adapted.fit(self.X_, self.y_) + return adapted + + def predict(self, X): + """ + Adapts the classifier to the target instances in `X`, and returns its label predictions for `X`. + Adaptation implies retraining. If you intend to generate inferences over many samples from the same + suspected distribution, then call the `get_classifier` method once, and use the returned classifier + over such samples. + + :param X: array-like of shape `(n_samples, n_features)` with the (unlabelled) target instances + :return: array-like of shape `(n_samples,)` with the class label predictions + """ + return self.get_classifier(X).predict(X) + + def predict_proba(self, X): + """ + Adapts the classifier to the target instances in `X`, and returns its posterior probabilities for `X`. + Adaptation implies retraining. If you intend to generate inferences over many samples from the same + suspected distribution, then call the `get_classifier` method once, and use the returned classifier + over such samples. + + :param X: array-like of shape `(n_samples, n_features)` with the (unlabelled) target instances + :return: array-like of shape `(n_samples, n_classes)` with posterior probabilities + """ + return self.get_classifier(X).predict_proba(X) diff --git a/quapy/method/_liep.py b/quapy/method/_liep.py new file mode 100644 index 0000000..3c28203 --- /dev/null +++ b/quapy/method/_liep.py @@ -0,0 +1,227 @@ +import warnings + +import numpy as np +from sklearn.base import BaseEstimator + +import quapy.functional as F +from quapy.method._helper import _labels_to_indices +from quapy.method.aggregative import AggregativeSoftQuantifier + + +def _normalize(v, eps=1e-12): + v = np.maximum(np.asarray(v, dtype=float), eps) + return v / v.sum() + + +def _bayes_prior_update(probs, new_prior, source_prior, eps=1e-12): + """ + Applies the standard label-shift prior correction: + + p_new(y|x) \\propto p_old(y|x) * new_prior(y) / source_prior(y) + + :param probs: array of shape `(n_instances, n_classes)` with the (uncorrected) posterior probabilities + :param new_prior: array of shape `(n_classes,)`, the current estimate of the target prior + :param source_prior: array of shape `(n_classes,)`, the source (training) prior + :param eps: numerical stabilizer + :return: array of shape `(n_instances, n_classes)` with the corrected posterior probabilities + """ + weights = _normalize(new_prior, eps) / np.maximum(_normalize(source_prior, eps), eps) + updated = probs * weights[np.newaxis, :] + updated /= np.maximum(updated.sum(axis=1, keepdims=True), eps) + return updated + + +def _confusion_statistic(y_true, y_pred, n_classes, mode='recall', eps=1e-12): + """ + Computes the minimum diagonal statistic of the confusion matrix, used to choose the fraction of target + instances to be retained in the high-confidence set (LEIP, Section 4.1). + + :param y_true: array of true label indices in `[0, n_classes)` + :param y_pred: array of predicted label indices in `[0, n_classes)` + :param n_classes: number of classes + :param mode: "recall" (per-true-class recall, i.e., diag(C) / row sums) or "precision" (per-predicted-class + precision, i.e., diag(C) / column sums) + :param eps: numerical stabilizer + :return: float, the minimum, across classes with at least one valid instance, of the chosen statistic + """ + C = np.zeros((n_classes, n_classes), dtype=float) + for yt, yp in zip(y_true, y_pred): + C[yt, yp] += 1. + + denom = C.sum(axis=1) if mode == 'recall' else C.sum(axis=0) + diag = np.diag(C) + valid = denom > 0 + if not np.any(valid): + raise ValueError('no valid classes found in the confusion matrix') + + return float(np.min(diag[valid] / np.maximum(denom[valid], eps))) + + +class LEIP(AggregativeSoftQuantifier): + """ + `Label Shift Estimation With Incremental Prior update + `_ (LEIP). + + Zhang, Y., Batista, G., & Kanhere, S.S. (2025). Label Shift Estimation With Incremental Prior Update. + In Proceedings of the 2025 SIAM International Conference on Data Mining (SDM), pp. 134-142. + + LEIP is a variant of :class:`EMQ` (aka SLD) that estimates the target prior without resorting to an + iterative Expectation-Maximization procedure. The target instances are first split into a high-confidence + set (those instances whose maximum posterior probability is at or above a threshold `tau`) and a + low-confidence set. The high-confidence set is used to obtain an initial classify-and-count estimate of + the target prior. The low-confidence instances are then relabelled one at a time, in order of decreasing + confidence, applying at each step the standard label-shift Bayes correction (see the module-level function + :func:`_bayes_prior_update`) using the running target-prior estimate; each newly assigned pseudo-label + updates this running estimate before the next instance is processed. Finally, a single Bayes correction is + applied to the whole target set using the resulting prior, and a hard-label classify-and-count pass on the + corrected posteriors yields the returned prevalence estimate. + + :param classifier: a scikit-learn's BaseEstimator, or None, in which case the classifier is taken to be + the one indicated in `qp.environ['DEFAULT_CLS']` + + :param fit_classifier: whether to train the classifier (default is True). Set to False if the + given classifier has already been trained. + + :param val_split: specifies the data used, when `tau=None`, for estimating the fraction of target + instances to retain in the high-confidence set (see `threshold_mode`). This specification can be made + as a float in (0, 1) indicating the proportion of stratified held-out validation set to be extracted + from the training set; or as an integer (default 5), indicating that the predictions are to be + generated in a `k`-fold cross-validation manner (with this integer indicating the value for `k`); or + as a tuple `(X, y)` defining the specific set of data to use for validation. This parameter is ignored + when a fixed `tau` is provided. + + :param tau: float or None (default). A fixed confidence threshold used to split the target instances into + the high- and low-confidence sets. If None, the threshold is instead derived, independently for each + target sample, as the quantile of that sample's own posterior-confidence distribution corresponding to + a retention fraction estimated from `val_split` (see `threshold_mode`). + + :param threshold_mode: either "recall" (default) or "precision", the statistic computed on the validation + confusion matrix that determines the fraction of target instances to be retained in the high-confidence + set (see the module-level function :func:`_confusion_statistic`). Only used when `tau=None`. + + :param count_smoothing: float, additive smoothing applied to the pseudo-label counts used to compute the + running target-prior estimate (default 0.0, matching the original paper; a small positive value, e.g., + 1e-8, can instead be used for extra numerical robustness). + """ + + THRESHOLD_MODES = ['recall', 'precision'] + + def __init__(self, classifier: BaseEstimator = None, fit_classifier=True, val_split=5, tau=None, + threshold_mode='recall', count_smoothing=0.0): + + assert threshold_mode in LEIP.THRESHOLD_MODES, \ + f'invalid {threshold_mode=}; valid ones are {LEIP.THRESHOLD_MODES}' + assert tau is None or isinstance(tau, (int, float)), f'invalid {tau=}; must be None or a number' + assert isinstance(count_smoothing, (int, float)) and count_smoothing >= 0, \ + f'invalid {count_smoothing=}; must be a non-negative number' + + super().__init__(classifier, fit_classifier, val_split) + self.tau = tau + self.threshold_mode = threshold_mode + self.count_smoothing = count_smoothing + + def _check_init_parameters(self): + if self.tau is None: + if self.val_split is None: + raise ValueError( + 'LEIP requires validation predictions to estimate the retention threshold when tau=None; ' + 'please set val_split to an integer, float, or validation tuple, or provide a fixed tau.' + ) + elif self.val_split is not None: + warnings.warn(f'a fixed {self.tau=} was set; {self.val_split=} will be ignored, since LEIP only ' + f'needs validation data to estimate a retention threshold when tau is not fixed.') + + def classify(self, X): + """ + Provides the posterior probabilities for the given instances. + + :param X: array-like of shape `(n_instances, n_dimensions,)` + :return: np.ndarray of shape `(n_instances, n_classes,)` with posterior probabilities + """ + return self.classifier.predict_proba(X) + + def classifier_fit_predict(self, X, y): + classif_predictions = super().classifier_fit_predict(X, y) + self.train_prevalence = F.prevalence_from_labels(y, classes=self.classes_) + return classif_predictions + + def aggregation_fit(self, classif_predictions, labels): + """ + Trains the aggregation function of LEIP. When a fixed `tau` was not provided at construction time, + this estimates, from the validation predictions, the fraction of target instances that should be + retained in the high-confidence set (see :func:`_confusion_statistic`). + + :param classif_predictions: array-like with the posterior probabilities of the validation instances + :param labels: array-like with the true labels associated to each classifier prediction + """ + if self.tau is None: + n_classes = len(self.classes_) + y_val_idx = _labels_to_indices(labels, self.classes_) + y_val_pred_idx = classif_predictions.argmax(axis=1) + self.retain_fraction_ = _confusion_statistic( + y_val_idx, y_val_pred_idx, n_classes, mode=self.threshold_mode + ) + else: + self.retain_fraction_ = None + + def aggregate(self, classif_posteriors): + prevalence, _ = self._leip(classif_posteriors) + return prevalence + + def predict_proba(self, instances): + """ + Returns the posterior probabilities updated by the final Bayes correction applied by LEIP. + + :param instances: np.ndarray of shape `(n_instances, n_dimensions)` + :return: np.ndarray of shape `(n_instances, n_classes)` + """ + classif_posteriors = self.classify(instances) + _, corrected_posteriors = self._leip(classif_posteriors) + return corrected_posteriors + + def _leip(self, target_probs): + n_classes = len(self.classes_) + source_prior = self.train_prevalence + + target_conf = target_probs.max(axis=1) + target_top = target_probs.argmax(axis=1) + + if self.tau is not None: + tau = self.tau + else: + retain_fraction = float(np.clip(self.retain_fraction_, 0., 1.)) + if retain_fraction <= 0: + tau = np.inf + elif retain_fraction >= 1: + tau = -np.inf + else: + tau = np.quantile(target_conf, 1. - retain_fraction) + self.tau_ = tau + + # high-confidence set A: an initial classify-and-count estimate of the target prior + A_mask = target_conf >= tau + counts = np.full(n_classes, self.count_smoothing, dtype=float) + if np.any(A_mask): + counts += np.bincount(target_top[A_mask], minlength=n_classes) + current_prior = counts / counts.sum() + else: + current_prior = source_prior.copy() + + # low-confidence set B, relabelled incrementally in order of decreasing confidence + B_indices = np.where(~A_mask)[0] + B_indices = B_indices[np.argsort(-target_conf[B_indices])] + + for idx in B_indices: + corrected = _bayes_prior_update(target_probs[idx:idx + 1], current_prior, source_prior)[0] + pseudo_label = int(np.argmax(corrected)) + counts[pseudo_label] += 1. + current_prior = counts / counts.sum() + + self.intermediate_prior_ = current_prior.copy() + + # final full-batch Bayes correction and hard-label classify-and-count pass + corrected_posteriors = _bayes_prior_update(target_probs, current_prior, source_prior) + final_labels = corrected_posteriors.argmax(axis=1) + prevalence = F.prevalence_from_labels(final_labels, classes=np.arange(n_classes)) + + return prevalence, corrected_posteriors diff --git a/quapy/tests/test_labelshift.py b/quapy/tests/test_labelshift.py new file mode 100644 index 0000000..adaae03 --- /dev/null +++ b/quapy/tests/test_labelshift.py @@ -0,0 +1,84 @@ +import unittest + +import numpy as np +from sklearn.linear_model import LogisticRegression +from sklearn.exceptions import NotFittedError +from sklearn.utils.validation import check_is_fitted + +from quapy.classification.labelshift import LabelShiftedClassifier +from quapy.method.aggregative import ACC, BBSEhard, BBSEsoft, ImportanceWeightQuantifier, RLLS +from quapy.functional import check_prevalence_vector +from quapy.tests._synthetic import make_dataset + + +class TestLabelShift(unittest.TestCase): + + dataset = make_dataset(n_train=200, n_test=100, n_classes=3, n_features=10, random_state=7, name='labelshift') + + def test_importance_weight_quantifiers(self): + Xtr, ytr = self.dataset.training.Xy + Xte = self.dataset.test.X + for cls in [BBSEhard, BBSEsoft, RLLS]: + if cls is RLLS: + try: + import cvxpy # noqa: F401 + except ImportError: + continue + + q = cls(LogisticRegression(max_iter=2000), val_split=3) + q.fit(Xtr, ytr) + + weights = q.get_importance_weights(Xte) + prevalence, weights2 = q.quantify_and_weigh(Xte) + prevalence_direct = q.predict(Xte) + + self.assertTrue(check_prevalence_vector(prevalence)) + np.testing.assert_allclose(weights, weights2) + np.testing.assert_allclose(prevalence, prevalence_direct) + # no leftover mutable per-call state (safety against races under concurrent calls) + self.assertFalse(hasattr(q, 'weights')) + self.assertFalse(hasattr(q, 'last_w_')) + + def test_label_shifted_classifier_generic_quantifier(self): + Xtr, ytr = self.dataset.training.Xy + Xte = self.dataset.test.X + + base_classifier = LogisticRegression(max_iter=2000) + q = ACC(LogisticRegression(max_iter=2000), val_split=3) + self.assertNotIsInstance(q, ImportanceWeightQuantifier) + + lsc = LabelShiftedClassifier(base_classifier, q) + lsc.fit(Xtr, ytr) + + weights = lsc.get_importance_weights(Xte) + self.assertEqual(weights.shape, (self.dataset.training.n_classes,)) + self.assertTrue(np.all(weights >= lsc.weight_epsilon)) + + preds = lsc.predict(Xte) + probs = lsc.predict_proba(Xte) + self.assertEqual(len(preds), len(Xte)) + self.assertEqual(probs.shape, (len(Xte), self.dataset.training.n_classes)) + np.testing.assert_allclose(probs.sum(axis=1), 1.0, rtol=1e-5) + + # the classifier instance passed at construction time must be left untouched + with self.assertRaises(NotFittedError): + check_is_fitted(base_classifier) + + def test_label_shifted_classifier_importance_weight_quantifier(self): + Xtr, ytr = self.dataset.training.Xy + Xte = self.dataset.test.X + + q = BBSEhard(LogisticRegression(max_iter=2000), val_split=3) + self.assertIsInstance(q, ImportanceWeightQuantifier) + + lsc = LabelShiftedClassifier(LogisticRegression(max_iter=2000), q) + lsc.fit(Xtr, ytr) + + weights = lsc.get_importance_weights(Xte) + adapted = lsc.get_classifier(Xte) + check_is_fitted(adapted) # should not raise + self.assertEqual(set(adapted.get_params()['class_weight'].keys()), set(lsc.classes_)) + + +if __name__ == '__main__': + unittest.main() From 790299b9bb9c4c66998c598477da3c44b681ba90 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Mon, 5 Oct 2026 13:52:29 +0200 Subject: [PATCH 7/9] liep draft removed --- quapy/method/_liep_draft.py | 261 ------------------------------------ quapy/method/aggregative.py | 7 +- 2 files changed, 2 insertions(+), 266 deletions(-) delete mode 100644 quapy/method/_liep_draft.py diff --git a/quapy/method/_liep_draft.py b/quapy/method/_liep_draft.py deleted file mode 100644 index 62c9bf7..0000000 --- a/quapy/method/_liep_draft.py +++ /dev/null @@ -1,261 +0,0 @@ -import numpy as np - - -def normalize(v, eps=1e-12): - v = np.asarray(v, dtype=float) - v = np.maximum(v, eps) - return v / v.sum() - - -def bayes_prior_update(probs, new_prior, source_prior, eps=1e-12): - """ - Applies the standard label-shift prior correction: - - p_new(y|x) ∝ p_old(y|x) * new_prior(y) / source_prior(y) - - Parameters - ---------- - probs : array, shape (n_samples, n_classes) or (n_classes,) - Probabilistic outputs of the classifier. - new_prior : array, shape (n_classes,) - Current estimate of the target prior. - source_prior : array, shape (n_classes,) - Source/training prior. - """ - probs = np.asarray(probs, dtype=float) - one_dim = probs.ndim == 1 - if one_dim: - probs = probs[None, :] - - new_prior = normalize(new_prior, eps=eps) - source_prior = normalize(source_prior, eps=eps) - - weights = new_prior / np.maximum(source_prior, eps) - updated = probs * weights[None, :] - updated = updated / np.maximum(updated.sum(axis=1, keepdims=True), eps) - - return updated[0] if one_dim else updated - - -def confusion_statistic(y_true, y_pred, n_classes, mode="recall", eps=1e-12): - """ - Computes the minimum diagonal statistic used to choose tau. - - mode="recall": - diag(C) / row sums, i.e. per-true-class recall. - - mode="precision": - diag(C) / column sums, i.e. per-predicted-class precision. - - The paper calls this quantity 'minimum recall', but also describes - a column-normalized confusion matrix, which would correspond more - closely to precision. I expose both options. - """ - y_true = np.asarray(y_true, dtype=int) - y_pred = np.asarray(y_pred, dtype=int) - - C = np.zeros((n_classes, n_classes), dtype=float) - for yt, yp in zip(y_true, y_pred): - C[yt, yp] += 1.0 - - if mode == "recall": - denom = C.sum(axis=1) - elif mode == "precision": - denom = C.sum(axis=0) - else: - raise ValueError("mode must be 'recall' or 'precision'.") - - diag = np.diag(C) - valid = denom > 0 - - if not np.any(valid): - raise ValueError("No valid classes found in the confusion matrix.") - - scores = diag[valid] / np.maximum(denom[valid], eps) - return float(np.min(scores)) - - -def choose_tau_from_validation( - target_probs, - val_probs, - y_val, - mode="recall", - retain_fraction=None, -): - """ - Chooses tau following the spirit of Section 4.1. - - The paper says tau is selected from the top n percentile of target - confidences, with n = min-recall * 100. To make the selected set A - contain approximately min-recall fraction of the target samples, we set: - - tau = quantile(max_probs, 1 - min_recall) - - so that about min_recall of the target points satisfy max_prob >= tau. - - If retain_fraction is provided, it overrides the validation-derived value. - """ - target_probs = np.asarray(target_probs, dtype=float) - val_probs = np.asarray(val_probs, dtype=float) - y_val = np.asarray(y_val, dtype=int) - - n_classes = target_probs.shape[1] - y_val_pred = val_probs.argmax(axis=1) - - if retain_fraction is None: - retain_fraction = confusion_statistic( - y_true=y_val, - y_pred=y_val_pred, - n_classes=n_classes, - mode=mode, - ) - - retain_fraction = float(np.clip(retain_fraction, 0.0, 1.0)) - - target_conf = target_probs.max(axis=1) - - if retain_fraction <= 0: - tau = np.inf - elif retain_fraction >= 1: - tau = -np.inf - else: - tau = np.quantile(target_conf, 1.0 - retain_fraction) - - return tau, retain_fraction - - -def leip( - target_probs, - source_prior, - tau=None, - val_probs=None, - y_val=None, - threshold_mode="recall", - count_smoothing=0.0, - eps=1e-12, - return_details=False, -): - """ - LEIP: Label shift Estimation with Incremental Prior update. - - Parameters - ---------- - target_probs : array, shape (n_target, n_classes) - Probabilistic classifier outputs on the target/test set. - source_prior : array, shape (n_classes,) - Source class prior p_s(y). - tau : float or None - Confidence threshold. If None, it is estimated using validation data. - val_probs : array, shape (n_val, n_classes), optional - Validation probabilistic outputs, required if tau is None. - y_val : array, shape (n_val,), optional - Validation labels, required if tau is None. - threshold_mode : {"recall", "precision"} - Statistic used to choose tau from validation data. - count_smoothing : float - Optional additive smoothing for pseudo-label counts. Set to 0.0 - for a closer implementation of the paper; use a small value such - as 1e-8 for extra numerical robustness. - eps : float - Numerical stabilizer. - return_details : bool - If True, returns diagnostic information. - - Returns - ------- - estimated_prior : array, shape (n_classes,) - Estimated target class distribution. - details : dict, optional - Returned only if return_details=True. - """ - target_probs = np.asarray(target_probs, dtype=float) - if target_probs.ndim != 2: - raise ValueError("target_probs must have shape (n_samples, n_classes).") - - n_target, n_classes = target_probs.shape - source_prior = normalize(source_prior, eps=eps) - - if n_classes != len(source_prior): - raise ValueError("source_prior must have one entry per class.") - - # Step 1: choose tau if needed - if tau is None: - if val_probs is None or y_val is None: - raise ValueError("val_probs and y_val are required when tau is None.") - - tau, retain_fraction = choose_tau_from_validation( - target_probs=target_probs, - val_probs=val_probs, - y_val=y_val, - mode=threshold_mode, - ) - else: - retain_fraction = None - - target_conf = target_probs.max(axis=1) - target_top = target_probs.argmax(axis=1) - - # Step 2: high-confidence set A - A_mask = target_conf >= tau - A_labels = target_top[A_mask] - - counts = np.full(n_classes, count_smoothing, dtype=float) - - if len(A_labels) > 0: - counts += np.bincount(A_labels, minlength=n_classes) - current_prior = counts / counts.sum() - else: - # Fallback if tau is too strict. - # One could also use classify-and-count over the full target set. - current_prior = source_prior.copy() - - # Step 3: low-confidence set B, sorted by decreasing confidence - B_indices = np.where(~A_mask)[0] - B_indices = B_indices[np.argsort(-target_conf[B_indices])] - - # Step 4: incremental pass over B - incremental_labels = [] - - for idx in B_indices: - corrected = bayes_prior_update( - probs=target_probs[idx], - new_prior=current_prior, - source_prior=source_prior, - eps=eps, - ) - pseudo_label = int(np.argmax(corrected)) - incremental_labels.append(pseudo_label) - - counts[pseudo_label] += 1.0 - current_prior = counts / counts.sum() - - estimated_intermediate_prior = current_prior.copy() - - # Step 5: final complete pass over all target instances - corrected_all = bayes_prior_update( - probs=target_probs, - new_prior=estimated_intermediate_prior, - source_prior=source_prior, - eps=eps, - ) - - final_labels = corrected_all.argmax(axis=1) - estimated_prior = np.bincount(final_labels, minlength=n_classes).astype(float) - estimated_prior /= estimated_prior.sum() - - if not return_details: - return estimated_prior - - details = { - "tau": tau, - "retain_fraction": retain_fraction, - "n_A": int(A_mask.sum()), - "n_B": int((~A_mask).sum()), - "A_mask": A_mask, - "intermediate_prior": estimated_intermediate_prior, - "final_labels": final_labels, - "corrected_probs": corrected_all, - "incremental_labels": np.asarray(incremental_labels, dtype=int), - } - - return estimated_prior, details \ No newline at end of file diff --git a/quapy/method/aggregative.py b/quapy/method/aggregative.py index 1decb39..6dcd025 100644 --- a/quapy/method/aggregative.py +++ b/quapy/method/aggregative.py @@ -2218,11 +2218,8 @@ KDEyML = _kdey.KDEyML KDEyHD = _kdey.KDEyHD KDEyCS = _kdey.KDEyCS -try: - from . import _liep_draft as _liep - LEIP = _liep.LEIP -except AttributeError: - LEIP = "LEIP is not available (incomplete implementation in _liep_draft.py)" +from . import _liep +LEIP = _liep.LEIP # --------------------------------------------------------------- From e4282502610337cf92f53bb47bb1ace4dd62bd99 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Mon, 5 Oct 2026 14:43:20 +0200 Subject: [PATCH 8/9] GMNet added --- CHANGE_LOG.txt | 7 ++- README.md | 5 +- examples/22.histnet_gmnet_bagtraining.py | 68 ++++++++++++++++++++++++ quapy/__init__.py | 2 +- quapy/method/_gmnet.py | 2 +- quapy/method/_histnet.py | 2 +- 6 files changed, 79 insertions(+), 7 deletions(-) create mode 100644 examples/22.histnet_gmnet_bagtraining.py diff --git a/CHANGE_LOG.txt b/CHANGE_LOG.txt index a364696..370a868 100644 --- a/CHANGE_LOG.txt +++ b/CHANGE_LOG.txt @@ -1,7 +1,12 @@ +Change Log 0.2.3 +----------------- + +- Added GMNet, based on the original implementation https://github.com/AICGijon/gmnet + Change Log 0.2.2 ----------------- -- Added HistNetQ, based on the original implementation https://github.com/pglez84/histnetq +- Added HistNetQ, based on the original implementation https://github.com/pglez82/histnetq - Minor fixes diff --git a/README.md b/README.md index 4110376..0c20e90 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ # QuaPy -## version 0.2.1 +## version 0.2.3 QuaPy is an open source framework for quantification (a.k.a. supervised prevalence estimation, or learning to quantify) written in Python. @@ -15,9 +15,8 @@ for facilitating the analysis and interpretation of the experimental results. ### Last updates: -* Version 0.2.1 is released! major changes can be consulted [here](CHANGE_LOG.txt). +* Version 0.2.3 is released! major changes can be consulted [here](CHANGE_LOG.txt). * The developer API documentation is available [here](https://hlt-isti.github.io/QuaPy/index.html) -* Manuals are available [here](https://hlt-isti.github.io/QuaPy/manuals.html) ### Installation diff --git a/examples/22.histnet_gmnet_bagtraining.py b/examples/22.histnet_gmnet_bagtraining.py new file mode 100644 index 0000000..985c138 --- /dev/null +++ b/examples/22.histnet_gmnet_bagtraining.py @@ -0,0 +1,68 @@ +import quapy as qp +from quapy.data.datasets import fetch_image_embeddings +from quapy.method.meta import HistNetQ, GMNet +from quapy.protocol import UPP + +# This example showcases HistNetQ and GMNet, two neural quantifiers that -- unlike the +# classify-then-aggregate methods seen in most other examples (e.g., EMQ, PACC) -- are trained +# end-to-end on whole samples ("bags") of known prevalence, instead of on individually labelled +# instances. Both require `torch` to be installed (GMNet additionally requires `geotorch`); see +# quapy.method._neural_bags.BagTrainedQuantifier, the base class they share, and the "Methods" manual +# for further details. +# +# In particular, we exercise `fit_from_samples`, the entry point meant for the situation in which +# individually labelled training instances are NOT available at all, and all one has access to is a +# pool of pre-built bags of known prevalence -- e.g., as released for the LeQua challenges, or as +# collected by a third party. We simulate such a pool here using CIFAR10 (see also example 20), by +# drawing a fixed set of bags from its (otherwise fully instance-labelled) training collection with +# the UPP protocol, and then training HistNetQ and GMNet as if that pool of bags -- and not the +# underlying instance-level labels -- were all we had access to. +# +# The settings below (number of bags, training epochs, patience) are deliberately modest so that the +# example runs in a reasonable time on a CPU; scale them up for serious experimentation. + +if __name__ == '__main__': + + qp.environ['SAMPLE_SIZE'] = 500 + + # The datasets available in quapy do not consist of raw image files, but are instead + # pre-generated embeddings (see the manuals for further information); 'features' corresponds to + # the penultimate-layer representations of a resnet18, which is what HistNetQ and GMNet expect to + # receive, since both operate on continuous instance representations. + print('fetching cifar10 embeddings') + train, test = fetch_image_embeddings(dataset_name='cifar10', embedding='features').train_test + print('training:', train) + print('test:', test) + + # Let us pretend that, instead of having access to the individually labelled instances of `train`, + # we are only given a pool of 1000 pre-built bags of 500 instances each, with heterogeneous class + # prevalence values (sampled uniformly at random from the simplex). This is exactly what + # fit_from_samples expects: an AbstractProtocol that already yields (sample, prevalence) pairs. + # UPP is deterministic by default (fixed random_state=0), so this pool is the same every time it + # is iterated, as befits a fixed, pre-built collection of bags. + print('simulating a pool of 1000 pre-built training bags (as if released by a third party)') + bag_pool = UPP(train, repeats=1000, random_state=0) + + # Neither method requires a classifier, nor (in this case) a feature extraction module: since the + # CIFAR10 embeddings are already vectorized, the default identity module is used, and both + # networks operate directly on the 512-dimensional resnet18 features + models = { + 'HistNetQ': HistNetQ(n_bins=8, bag_size=qp.environ['SAMPLE_SIZE'], train_epochs=50, patience=10, verbose=1), + 'GMNet': GMNet(n_gm_layers=8, num_gaussians=10, gaussian_dimensions=5, hidden_size_fe=(50,), + bag_size=qp.environ['SAMPLE_SIZE'], train_epochs=50, patience=10, verbose=1), + } + + for name, model in models.items(): + print(f'\nfitting {name} from the pre-built bag pool via fit_from_samples ' + f'(no instance-level labels are used)') + model.fit_from_samples(bag_pool) + + # we now evaluate both quantifiers on cifar10's test set, using a fresh artificial-prevalence + # protocol (note that, at prediction time, HistNetQ and GMNet behave just like any other + # quantifier: the bag-based training is an internal detail of these two particular methods) + test_prot = UPP(test, repeats=1000) + + for name, model in models.items(): + report = qp.evaluation.evaluation_report(model, protocol=test_prot, error_metrics=['mae']) + print(f'\n{name}:') + print(report.mean(numeric_only=True)) diff --git a/quapy/__init__.py b/quapy/__init__.py index 4ed6726..de4575c 100644 --- a/quapy/__init__.py +++ b/quapy/__init__.py @@ -17,7 +17,7 @@ try: except ImportError: plot = None -__version__ = '0.2.2' +__version__ = '0.2.3' def _default_cls(): diff --git a/quapy/method/_gmnet.py b/quapy/method/_gmnet.py index 195681d..6222705 100644 --- a/quapy/method/_gmnet.py +++ b/quapy/method/_gmnet.py @@ -1,7 +1,7 @@ """ GMNet implementation. -Ported from the reference implementation at https://github.com/pglez84/gmnet (the `GMNet`/ +Ported from the reference implementation at https://github.com/AICGijon/gmnet (the `GMNet`/ `DLQuantification` classes in that repo), adapted to QuaPy's own protocol-based sample generation (replacing that repo's custom, `quantificationlib`-backed bag generators), and reusing the shared bag-based training loop already factored out for :class:`quapy.method._histnet.HistNetQ` (see diff --git a/quapy/method/_histnet.py b/quapy/method/_histnet.py index ad58725..d718464 100644 --- a/quapy/method/_histnet.py +++ b/quapy/method/_histnet.py @@ -98,7 +98,7 @@ class _SigmoidHistogram(nn.Module): class HistNetQ(BagTrainedQuantifier): """ - Implementation of `HistNetQ `_, a neural network for + Implementation of `HistNetQ `_, a neural network for quantification that learns a differentiable histogram-based representation of a sample, trained end-to-end by minimizing a quantification loss over many samples ("bags") of known prevalence. The method was proposed in `Pérez-Mon, O., Moreo, A., Coz, JJ del, & González, P. (2025). From 5d3c2bef16b47272dbd1d09c14ea0f7e7a2daf38 Mon Sep 17 00:00:00 2001 From: Alejandro Moreo Date: Mon, 5 Oct 2026 14:46:51 +0200 Subject: [PATCH 9/9] histnet and gmnet manuals update --- docs/source/manuals/methods.md | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/docs/source/manuals/methods.md b/docs/source/manuals/methods.md index 26fea42..dcab373 100644 --- a/docs/source/manuals/methods.md +++ b/docs/source/manuals/methods.md @@ -907,7 +907,8 @@ estim_prevalence = model.predict(dataset.test.X) ### HistNetQ -QuaPy offers an implementation of HistNetQ, a deep learning model based on a differentiable +QuaPy offers an implementation of HistNetQ (based on the [original implementation](https://github.com/pglez82/histnetq)), +a deep learning model based on a differentiable histogram representation, presented in: [_Pérez-Mon, O., Moreo, A., del Coz, J.J., & González, P. (2025). @@ -936,7 +937,8 @@ method; see the API documentation for further details. ### GMNet -QuaPy offers an implementation of GMNet, a deep learning model that represents each instance of a +QuaPy offers an implementation of GMNet (based on the [original implementation](https://github.com/AICGijon/gmnet)), +a deep learning model that represents each instance of a bag by its likelihood under one or more learned mixtures of Gaussians, presented in: [_Pérez-Mon, O., del Coz, J.J., & González, P. (2026).