Exact R2 Indicator

README.md

@@ -15,10 +15,11 @@
     [Manuel López-Ibáñez](https://lopez-ibanez.eu),
     [Carlos M. Fonseca](https://eden.dei.uc.pt/~cmfonsec/),
     [Luís Paquete](https://eden.dei.uc.pt/~paquete/),
-    Andreia P. Guerreiro
-    Mickaël Binois.
+    Andreia P. Guerreiro,
+    Mickaël Binois,
     Leonardo C.T. Bezerra,
-    Fergus Rooney.
+    Fergus Rooney,
+    [Lennart Schäpermeier](https://schaepermeier.github.io).
 
 ---------------------------------------
 

bibkeys.txt

@@ -17,6 +17,7 @@ GruFon2009:emaa
 Grunert01
 GueFon2017hv4d
 GueFonPaq2021hv
+HanJas1998
 HerWer1987tabucol
 IshMasTanNoj2015igd
 Jen03
@@ -31,6 +32,7 @@ Molchanov2005theory
 Muller1959sphere
 RubMel1998simulation
 SchEsqLarCoe2012tec
+SchKer2025r2v2
 VelLam1998gp
 WuAza2001metrics
 ZhoZhaJin2009igdx

c/Makefile

@@ -62,6 +62,7 @@ SRCS =  avl.c                                                                \
         eafdiff.c                                                            \
         eaf_main.c                                                           \
         epsilon.c                                                            \
+        r2_exact.c                                                           \
         hv3dplus.c                                                           \
         hv4d.c                                                               \
         hv.c                                                                 \
@@ -91,6 +92,7 @@ HDRS =  avl.h                                                                \
         cvector.h                                                            \
         eaf.h                                                                \
         epsilon.h                                                            \
+        r2_exact.h                                                           \
         gcc_attribs.h                                                        \
         hv.h                                                                 \
         hvapprox.h                                                           \

c/r2_exact.c

@@ -0,0 +1,94 @@
+#include <stdlib.h>
+#include <limits.h>
+#include <float.h>
+#include "common.h"
+#include "sort.h"
+
+// Computes (two times) the utility of axis-parallel segment between y1, y2, and y2p.
+static double _utility(const double y1, const double y2, const double y2p)
+{
+   if (y1 == 0) {
+        return 0;
+    }
+
+    double w = y2 / (y1 + y2);
+    double wp;
+    if (y2p == DBL_MAX) {
+        wp = 1;
+    } else {
+        wp = y2p / (y1 + y2p);
+    }
+
+    return y1 * (wp * wp - w * w);
+}
+
+double r2_exact(const double * restrict data, size_t n, dimension_t dim,
+                const double * restrict ref)
+{
+    if (unlikely(n == 0)) return -1;
+
+    ASSUME(dim == 2);
+
+    // p is sorted by f1 (primarily), then f2 (secondarily)
+    const double **p = generate_sorted_doublep_2d(data, &n, DBL_MAX);
+
+    if (unlikely(!p)) return -1;
+
+    size_t j = 0;
+
+    // skip points not dominated by "ideal" ref point
+    while (j < n && p[j][0] < ref[0]) {
+        j++;
+    }
+
+    // no points to evaluate
+    if (unlikely(j == n)) {
+        if (p[j - 1][1] <= ref[1]) return 0.0; // ideal ref. is dominated.
+        return DBL_MAX; // ideal ref. is nondominated --> worst possible value
+    }
+
+    double prev_y1 = p[j][0] - ref[0];
+    double prev_y2 = p[j][1] - ref[1];
+
+    if (prev_y2 < 0) {
+        // ideal ref. is dominated.
+        return 0.0;
+    }
+
+    // first element
+    double r2_exact = _utility(prev_y1, prev_y2, DBL_MAX);
+    // printf("y2 segment: %f, %f, MAX: %f\n", p[j][0] - ref[0], p[j][1] - ref[1], _utility(p[j][0] - ref[0], p[j][1] - ref[1], DBL_MAX));
+
+    while(j < n - 1) {
+        j++;
+        double y1 = p[j][0] - ref[0];
+        double y2 = p[j][1] - ref[1];
+
+        // skip anything that's not dominated by ref
+        if (y2 < 0) continue;
+
+        // printf("%f\n", y1);
+        // printf("%f\n", y2);
+
+        if (y2 < prev_y2) {
+            // pj not dominated by previous non-dominated pj
+            r2_exact += _utility(prev_y2, prev_y1, y1) + _utility(y1, y2, prev_y2);
+
+            // printf("y1 segment: %f, %f, %f: %f\n", prev_y2, prev_y1, y1, _utility(prev_y2, prev_y1, y1));
+            // printf("y2 segment: %f, %f, %f: %f\n", y1, y2, prev_y2, _utility(y1, y2, prev_y2));
+
+            prev_y1 = y1;
+            prev_y2 = y2;
+        }
+    }
+
+    // last element
+    r2_exact += _utility(prev_y2, prev_y1, DBL_MAX);
+    // printf("y1 segment: %f, %f, MAX: %f\n", prev_y2, prev_y1, _utility(prev_y2, prev_y1, DBL_MAX));
+
+    // we omitted a 1/2 in the computation thus far:
+    r2_exact = 0.5 * r2_exact;
+
+    free(p);
+    return r2_exact;
+}

c/r2_exact.h

@@ -0,0 +1,31 @@
+/*************************************************************************
+
+ r2_exact.h
+
+ ---------------------------------------------------------------------
+
+                          Copyright (C) 2026
+        Lennart Schäpermeier <lennart.schaepermeier@uni-muenster.de>
+
+ This Source Code Form is subject to the terms of the Mozilla Public
+ License, v. 2.0. If a copy of the MPL was not distributed with this
+ file, You can obtain one at https://mozilla.org/MPL/2.0/.
+
+ ----------------------------------------------------------------------
+
+
+*************************************************************************/
+#ifndef R2_EXACT_H_
+#define R2_EXACT_H_
+
+#include <stdbool.h>
+#include "libmoocore-config.h"
+
+// C++ needs to know that types and declarations are C, not C++.
+BEGIN_C_DECLS
+
+MOOCORE_API double r2_exact(const double * restrict data, size_t n, dimension_t d, const double * restrict ref);
+
+END_C_DECLS
+
+#endif // R2_EXACT_H_

python/README.md

@@ -14,7 +14,8 @@
 
 **Contributors:**
     [Manuel López-Ibáñez](https://lopez-ibanez.eu),
-    Fergus Rooney.
+    Fergus Rooney,
+    [Lennart Schäpermeier](https://schaepermeier.github.io).
 
 ---------------------------------------
 

python/doc/source/REFERENCES.bib

@@ -361,6 +361,17 @@ @article{SchEsqLarCoe2012tec
   doi = {10.1109/TEVC.2011.2161872}
 }
 
+@article{SchKer2025r2v2,
+  author = {Sch{\"a}permeier, Lennart and  Pascal Kerschke },
+  title = {{R2} v2: The {Pareto}-compliant {R2} Indicator for Better
+                  Benchmarking in Bi-objective Optimization},
+  journal = {Evolutionary Computation},
+  year = 2025,
+  pages = {1--17},
+  month = oct,
+  doi = {10.1162/evco.a.366}
+}
+
 @article{WuAza2001metrics,
   author = {J. Wu and S. Azam},
   title = {Metrics for Quality Assessment of a Multiobjective Design
@@ -604,6 +615,18 @@ @incollection{Grunert01
                   paper.}
 }
 
+@techreport{HanJas1998,
+  author = { Michael Pilegaard Hansen  and  Andrzej Jaszkiewicz },
+  title = {Evaluating the quality of approximations to the non-dominated
+                  set},
+  institution = {Institute of Mathematical Modelling, Technical University of
+                  Denmark},
+  year = 1998,
+  number = {IMM-REP-1998-7},
+  address = {Lyngby, Denmark},
+  annote = {Proposed R2 indicator}
+}
+
 @incollection{IshMasTanNoj2015igd,
   booktitle = { Evolutionary Multi-criterion Optimization, EMO 2015 Part {I}},
   address = { Heidelberg, Germany},

python/doc/source/index.rst

@@ -28,7 +28,7 @@ The goal of the **moocore** project (https://github.com/multi-objective/moocore/
 
 * :ref:`Generate and transform nondominated sets <read_generate_transform_sets>`.
 * :ref:`Identify and filter dominated vectors <identifying_and_filtering_dominated_vectors>`.
-* :ref:`Quality metrics <unary_quality_metrics>` such as (weighted) hypervolume, epsilon, IGD+, etc.
+* :ref:`Quality metrics <unary_quality_metrics>` such as (weighted) hypervolume, epsilon, IGD+, (exact) R2, etc.
 * :ref:`Computation of the Empirical Attainment Function <empirical_attainment_function>`. The empirical attainment function (EAF) describes the probabilistic distribution of the outcomes obtained by a stochastic algorithm in the objective space.
 
 Most critical functionality is implemented in C, with the R and Python packages providing convenient interfaces to the C code.

python/doc/source/reference/functions.metrics.rst

@@ -190,6 +190,46 @@ Approximating the hypervolume metric
 Computing the hypervolume can be time consuming, thus several approaches have been proposed in the literature to approximate its value via Monte-Carlo or quasi-Monte-Carlo sampling :cite:p:`DenZha2019approxhv`. These methods are implemented in :func:`whv_hype` and :func:`hv_approx`.
 
 
+.. _r2_indicator:
+
+(Exact) R2 Indicator
+=============================================================================
+
+.. autosummary::
+   :toctree: generated/
+
+   r2_exact
+
+The unary R2 indicator is a quality indicator for a set :math:`A \subset \mathbb{R}^m`
+w.r.t. an ideal or utopian reference point :math:`\vec{r} \in \mathbb{R}^m`.
+It was originally proposed by :cite:t:`HanJas1998` and is defined as the
+expected Tchebycheff utility under a uniform distribution of weight vectors
+(w.l.o.g. assuming minimization):
+
+.. math::
+   R2(A) := \int_{w \in W} \min_{a \in A} \left\{ \max_{i=1,\dots,m} w_i (a_i - r_i) \right\} \, dw,
+
+where :math:`W` denotes the uniform distribution across weights:
+
+.. math::
+   W = \{w \in \mathbb R^m \mid w_i \geq 0, \sum_{i=1}^m w_i = 1\}.
+
+The R2 indicator is to be minimized and has an optimal value of 0 when
+:math:`\vec{r} \in A`.
+
+The exact R2 indicator is strongly Pareto-compliant, i.e., compatible with Pareto-optimality:
+
+.. math::
+   \forall A, B \subset \mathbb R^m: A \prec B \Rightarrow R2(A) < R2(B).
+
+Given an ideal or utopian reference ponint, which is available in most scenarios,
+all non-dominated solutions always contribute to the value of the exact R2 indicator.
+However, it is scale-dependent and care should be taken such that all objectives
+contribute approximately equally to the indicator, e.g., by normalizing the
+Pareto front to the unit hypercube.
+
+:func:`r2_exact` computes the exact R2 indicator for bi-objective solution
+sets in :math:`\mathcal O(N \log N)` following :cite:t:`SchKer2025r2v2`.
 
 Bibliography
 ============

python/doc/source/whatsnew/index.rst

@@ -8,7 +8,7 @@ Version 0.3.0 (development)
 ---------------------------
 
 - `moarchiving <https://cma-es.github.io/moarchiving/>`_ is now included in the benchmarks.
-
+- :func:`~moocore.r2_exact` implements the exact computation of the R2 indicator for bi-objective solution sets.
 
 Version 0.2.0 (10/01/2026)
 ---------------------------

python/examples/plot_metrics.py

@@ -98,10 +98,11 @@
 # %%
 #
 # However, both :func:`moocore.igd` and :func:`moocore.avg_hausdorff_dist`
-# incorrectly measure B as better than A, whereas :func:`moocore.igd_plus` and
-# :func:`moocore.hypervolume` correctly measure A as better than B (remember
-# that hypervolume must be maximized) and :func:`moocore.epsilon_additive` measures
-# both as equally good (epsilon is weakly Pareto compliant).
+# incorrectly measure B as better than A, whereas :func:`moocore.igd_plus`,
+# :func:`moocore.r2_exact` and :func:`moocore.hypervolume` correctly measure A
+# as better than B (remember that hypervolume must be maximized) and
+# :func:`moocore.epsilon_additive` measures both as equally good
+# (epsilon is weakly Pareto compliant).
 #
 
 pd.DataFrame(
@@ -112,14 +113,16 @@
             moocore.igd_plus(A, ref),
             moocore.epsilon_additive(A, ref),
             moocore.hypervolume(A, ref=10),
+            moocore.r2_exact(A, ref=0),
         ],
         B=[
             moocore.igd(B, ref),
             moocore.avg_hausdorff_dist(B, ref),
             moocore.igd_plus(B, ref),
             moocore.epsilon_additive(B, ref),
             moocore.hypervolume(B, ref=10),
+            moocore.r2_exact(B, ref=0),
         ],
     ),
-    index=["IGD", "Hausdorff", "IGD+", "eps+", "HV"],
+    index=["IGD", "Hausdorff", "IGD+", "eps+", "HV", "Exact R2"],
 )

python/pyproject.toml

@@ -23,6 +23,7 @@ license-files = [ "COPYRIGHTS" ]
 authors = [
   { name = "Manuel López-Ibáñez", email = "manuel.lopez-ibanez@manchester.ac.uk" },
   { name = "Fergus Rooney", email = "fergus.rooney@outlook.com" },
+  { name = "Lennart Schäpermeier", email = "lennart.schaepermeier@uni-muenster.de" },
 ]
 requires-python = ">=3.10"
 classifiers = [

python/src/moocore/__init__.py

@@ -10,6 +10,7 @@
     eafdiff,
     epsilon_additive,
     epsilon_mult,
+    r2_exact,
     filter_dominated,
     filter_dominated_within_sets,
     generate_ndset,
@@ -69,6 +70,7 @@
     "largest_eafdiff",
     "normalise",
     "pareto_rank",
+    "r2_exact",
     "read_datasets",
     "total_whv_rect",
     "vorob_dev",

python/src/moocore/_ffi_build.py

@@ -24,6 +24,7 @@
 #include "nondominated.h"
 #include "epsilon.h"
 #include "eaf.h"
+#include "r2_exact.h"
 #include "whv.h"
 #include "whv_hype.h"
 #include "hvapprox.h"
@@ -33,6 +34,7 @@
     "eaf.c",
     "eaf3d.c",
     "eafdiff.c",
+    "r2_exact.c",
     "hv.c",
     "hvapprox.c",
     "hv3dplus.c",