Identifying Consistent Statements about Numerical Data with Dispersion-Corrected Subgroup Discovery
Mario Boley Bryan R. Goldsmith Luca M. Ghiringhelli Jilles Vreeken


3D Printing of Polymer Bonded Rare-Earth Magnets With a Variable Magnetic Compound Density for a Predefined Stray Field
    D. Suess
    December 30, 2023
=======================================================================================================================




Existing algorithms for subgroup discovery with numerical targets do not optimize the error or target variable dispersion of the groups they find. This often leads to unreliable or inconsistent statements about the data, rendering practical applications, especially in scientific domains, futile. Therefore, we here extend the optimistic estimator framework for optimal subgroup discovery to a new class of objective functions: we show how tight estimators can be computed efficiently for all functions that are determined by subgroup size (non-decreasing dependence), the subgroup median value, and a dispersion measure around the median (non-increasing dependence). In the important special case when dispersion is measured using the mean absolute deviation from the median, this novel approach yields a linear time algorithm. Empirical evaluation on a wide range of datasets shows that, when used within branch-and-bound search, this approach is highly efficient and indeed discovers subgroups with much smaller errors.




§ INTRODUCTION

Subgroup discovery is a well-established KDD technique (<cit.>; see <cit.> for a recent survey) with applications, e.g., in Medicine <cit.>, Social Science <cit.>, and Materials Science <cit.>. In contrast to global modeling, which is concerned with the complete characterization of some variable defined for a given population, subgroup discovery aims to detect intuitive descriptions or selectors of subpopulations in which, locally, the target variable takes on a useful distribution.
In scientific domains, like the ones mentioned above, 
such local patterns are typically considered useful if they are not too specific (in terms of subpopulation size) and indicate insightful facts about the underlying physical process that governs the target variable.
Such facts could for instance be: `patients of specific demographics experience a low response to some treatment' or `materials with specific atomic composition exhibit a high thermal conductivity'.
For numeric (metric) variables, subgroups need to satisfy two criteria to truthfully represent such statements: the local distribution of the target variable must have a shifted central tendency (effect), and group members must be described well by that shift (consistency). The second requirement is captured by the group's dispersion, which determines the average error of associating group members with the central tendency value <cit.>.


Despite all three parameters—size, central tendency, and dispersion—being important, the only known approach for the efficient discovery of globally optimal subgroups, branch-and-bound search <cit.>, is restricted to objective functions that only take into account size and central tendency.
That is, if we denote by Q some subpopulation of our global population P then the objective functions f currently available to branch-and-bound can be written as 

    f(Q)=g(Q, Q)

where  is some measure of central tendency (usually mean or median) and g is a function that is monotonically increasing in the subpopulation size Q.
A problem with all such functions is that they inherently favor larger groups with scattered target values over smaller more focused groups with the same central tendency. 
That is, they favor the discovery of inconsistent statements over consistent ones—surprisingly often identifying groups with a local error that is almost as high or even higher than the global error (see Fig. <ref> for an illustration of this problem that abounded from the authors' research in Materials Science).
Although dispersion-corrected objective functions that counter-balance size by dispersion have been proposed (e.g., `t-score' by  or `mmad' by ),
it remained unclear how to employ such functions outside of heuristic optimization frameworks such as greedy beam search <cit.> or selector sampling <cit.>.
Despite often finding interesting groups, such frameworks do not guarantee the detection of optimal results, which can not only be problematic for missing important discoveries but also because they therefore can never guarantee the absence of high quality groups—which often is an insight equally important as the presence of a strong pattern.
For instance, in our example in Fig. <ref>, it would be remarkable to establish that long-range interactions are to a large degree independent of nanocluster geometry.

Therefore, in this paper (Sec. <ref>), we extend branch-and-bound search to objective functions of the form

    f(Q)=g(Q, Q, Q)

where g is monotonically increasing in the subpopulation size, monotonically decreasing in any dispersion measure  around the median, and, besides that, depends only (but in arbitrary form) on the subpopulation median.
This involves developing an efficient algorithm for computing the tight optimistic estimator given by the optimal value of the objective function among all possible subsets of target values:

    f(Q)=max{f(R) R ⊆ Q} ,

which has been shown to be a crucial ingredient for the practical applicability of branch-and-bound <cit.>.
So far, the most general approach to this problem (first codified in <cit.>; generalized here in Sec. <ref>) is to maintain a sorted list of target values throughout the search process and then to compute Eq. (<ref>) as the maximum of all subsets R_i ⊆ Q that contain all target values of Q down to target value i—an algorithm that does not generalize to objective functions depending on dispersion.
This paper presents an alternative idea (Sec. <ref>)
where we do not fix the size of subset R_i as in the previous approach but instead fix its median to target value i. It turns out that this suffices to efficiently compute the tight optimistic estimator for all objective functions of the form of Eq. (<ref>).
Moreover, we end up with a linear time algorithm (Sec. <ref>) in the important special case where the dependence on size and dispersion is determined 
by the dispersion-corrected coverage defined by

    Q=Q/Pmax{1-Q/P,0}

where  denotes the mean absolute deviation from the median.
This is the same computational complexity as the objective function itself.
Consequently, this new approach can discover subgroups according to a more refined selection criterion without increasing the worst-case computational cost.
Additionally, as demonstrated by empirical results on a wide range of datasets (Sec. <ref>), it is also highly efficient and
successfully reduces the error of result subgroups in practice.



§ SUBGROUP DISCOVERY


Before developing the novel approach to tight optimistic estimator computation, we recall in this section the necessary basics of optimal subgroup discovery with numeric target attributes. We focus on concepts that are essential from the optimization point of view (see, e.g.,  and references therein for statistical considerations). As notional convention, we are using the symbol m for a positive integer m to denote the set of integers {1,…,m}.
Also, for a real-valued expression x we write (x)_+ to denote max{x,0}.
A summary of the most important notations used in this paper can be found in Appendix <ref>.



 §.§ Description languages, objective functions, and closed selectors


Let P denote our given global population of entities, for each of which we know the value of a real target variable y:P →ℝ and additional descriptive information that is captured in some abstract description language  of subgroup selectors σP{,}.
Each of these selectors describes a subpopulation σ⊆ P defined by 

    σ={p ∈ P σ(p)=}

that is referred to as the extension of σ.
Subgroup discovery is concerned with finding descriptions σ∈ that have a useful (or interesting) distribution of target values in their extension y_σ={y(p) : p ∈σ}.
This notion of usefulness is given by an objective function f. That is, the formal goal is to find elements σ∈ with maximal f(σ).
Since we assume f to be a function of the multiset of y-values, let us define f(σ)=f(σ)=f(y_σ) to be used interchangeably for convenience. 
One example of a commonly used objective function is the impact measure  (see ; here a scaled but order-equivalent version is given) defined by

    Q=Q ( Q-P/P-P )_+

where Q=Q/P denotes the coverage or relative size of Q (here—and wherever else convenient—we identify a subpopulation Q ⊆ P with the multiset of its target values).

The standard description language in the subgroup discovery literature[In this article we remain with this basic setting for the sake of simplicity. It is, however, useful to note that several generalizations of this concept have been proposed <cit.>, to which the contributions of this paper remain applicable.] is the language  consisting of logical conjunctions of a number of base propositions (or predicates).
That is, σ∈ are of the form

    σ(·) ≡π_i_1(·) ∧…∧π_i_l(·)

where the π_i_j are taken from a pool of base propositions Π={π_1,…,π_k}. These propositions usually correspond to equality or inequality constraints with respect to one variable x out of a set of description variables {x_1,…,x_n} that are observed for all population members (e.g., π(p)≡ x(p) ≥ v). However, for the scope of this paper it is sufficient to simply regard them as abstract Boolean functions πP{,}.
In this paper, we focus in particular on the refined language of closed conjunctions ⊆ <cit.>, which is 
defined as ={σ∈: σ=σ} by the fixpoints of the closure operation  given by 

    σ = ⋀{π∈Π: π⊇σ} .

These are selectors to which no further proposition can be added without reducing their extension, and it can be shown that  contains at most one selector for each possible extension.
While this can reduce the search space for finding optimal subgroups by several orders of magnitude, closed conjunctions are the longest (and most redundant) description for their extension and thus do not constitute intuitive descriptions by themselves. 
Hence, for reporting concrete selectors (as in Fig. <ref>), closed conjunctions have to be simplified to selectors of approximately minimum length that describe the same extension <cit.>.



 §.§ Branch-and-bound and optimistic estimators


bstbbBst-BB
qemptyempty
qaddalladdAll
qtoptop
sucsuc
argmaxargmax

The standard algorithmic approach for finding optimal subgroups with respect to a given objective function is branch-and-bound search—a versatile algorithmic puzzle solving framework with several forms and flavors <cit.>. At its core, all of its variants assume the availability and efficient computability of two ingredients: 


  * A refinement operator 2^ that is monotone, i.e., for σ, φ∈ with φ∈(σ) it holds that φ⊆σ,
and that non-redundantly generates .
That is, there is a root selector ∈ such that for every σ∈ there is a unique sequence of selectors =σ_0, σ_1, …, σ_l=σ with σ_i∈(σ_i-1). In other words, the refinement operator implicitly represents a directed tree (arborescence) on the description language  rooted in . 

  * An optimistic estimator (or bounding function) f that bounds from above the attainable subgroup value 
of a selector among all more specific selectors, i.e., it holds that f(σ) ≥ f(φ) for all φ∈ with φ⊆σ.

Based on these ingredients, a branch-and-bound algorithm simply enumerates all elements of  starting from  using  (branch), but—based on f—avoids expanding descriptions that cannot yield an improvement over the best subgroups found so far (bound). Depending on the order in which language elements are expanded, one distinguishes between depth-first, breadth-first, breadth-first iterating deepening, and best-first search.
In the last variant, the optimistic estimator is not only used for pruning the search space, but also to select the next element to be expanded, which is particularly appealing for informed, i.e., tight, optimistic estimators.
An important feature of branch-and-bound is that it effortlessly allows to speed-up the search in a sound way by relaxing the result requirement from being f-optimal to just being an a-approximation. That is, the found solution σ satisfies for all σ' ∈ that f(σ)/f(σ') ≥ a for some approximation factor a ∈ (0,1].
The pseudo-code given in Alg. <ref> summarizes all of the above ideas.
Note that, for the sake of clarity, we omitted here some other common parameters such as a depth-limit and multiple solutions (top-k), which are straightforward to incorporate <cit.>.

An efficiently computable refinement operator has to be constructed specifically for the desired description language. For example for the language of conjunctions , one can define _cnj by

    _cnj(σ)={σ∧π_i max{jπ_j ∈σ} < i ≤ k}

where we identify a conjunction with the set of base propositions it contains. 
For the closed conjunctions , 
let us define the lexicographical prefix of a conjunction σ∈ and a base proposition index i ∈k as σi=σ∩{π_1, …, π_i}. Moreover, let us denote with σ the minimal index such that the i-prefix of σ is extension-preserving, i.e.,
σ=min{i σi=σ}.
With this we can construct a refinement operator <cit.> _ccj as

    _ccj(σ)={φφ=σ∧π_j, σ < j ≤ k, π_j ∉σ, φj=σj} .

That is, a selector φ is among the refinements of σ if φ can be generated by an application of the closure operator given in Eq. (<ref>) that is prefix-preserving.

How to obtain an optimistic estimator for an objective function of interest depends on the definition of that objective.
For instance, the coverage function  is a valid optimistic estimator for the impact function  as defined in Eq. (<ref>), because the second factor of the impact function is upper bounded by 1.
In fact there are many different optimistic estimators for a given objective function.
Clearly, the smaller the value of the bounding function for a candidate subpopulation, the higher is the potential for pruning the corresponding branch from the enumeration tree.
Ideally, one would like to use f(σ)=max{f(φ) φ⊆σ}, which is the most strict function that still is a valid optimistic estimator.
In general, however, computing this function is as hard as the original subgroup optimization problem we started with.
Therefore, as a next best option, one can disregard selectability and consider the (selection-unaware) tight optimistic estimator given by

    f(σ)=max{f(R) R ⊆σ} .

This leaves us with a new combinatorial optimization problem: given a subpopulation Q ⊆ P, find a sub-selection of Q that maximizes f. In the following section we will discuss strategies for solving this optimization problem efficiently for different classes of objective functions—including dispersion-corrected objectives. 





§ EFFICIENTLY COMPUTABLE TIGHT OPTIMISTIC ESTIMATORS


We are going to develop an efficient algorithm for the tight optimistic estimator in three steps: First, we review and reformulate a general algorithm for the classic case of non-dispersion-aware objective functions. Then we transfer the main idea of this algorithm to the case of dispersion-corrected objectives based on the median, and finally we consider a subclass of these functions where the approach can be computed in linear time.
Throughout this section we will identify a given subpopulation Q ⊆ P with the multiset of its target values {y_1,…,y_m} and assume that the target values are indexed in ascending order, i.e., y_i ≤ y_j for i ≤ j.
Also, it is helpful to define the following partial order defined on finite multisets.
Let Y={y_1,…,y_m} and Z={z_1,…,z_m'} be two multisets such that their elements are indexed in ascending order. We say that Y is element-wise less or equal to Z and write Y  Z if y_i ≤ z_i for all i ∈min{m,m'}.



 §.§ The standard case: monotone functions of a central tendency measure


The most general previous approach for computing the tight optimistic estimator for subgroup discovery with a metric target variable is described by 
<cit.> where it is referred to as estimation by ordering. Here, we review this approach and give a uniform and generalized version of that paper's results.
For this, we define the general notion of a measure of central tendency as follows.

We call a mapping ^ a (monotone) measure of central tendency if for all multisets Y,Z ∈^ with Y  Z it holds that Y≤Z.

One can check that this definition applies to the standard measures of central tendency, i.e., the arithmetic and geometric mean as well as the median[In this paper, we are using the simple definition of the median as the 0.5-quantile (as opposed to defining it as (y_m/2+y_1+m/2)/2 for even m), which simplifies many of the definitions below and additionally is well-defined in settings where averaging of target values is undesired.] Q=y_⌈ m/2 ⌉, and also to weighted variants of them (note, however, that it does not apply to the mode).
With this we can define the class of objective functions for which the tight optimistic estimator can be computed efficiently by the standard approach as follows.
We call fP a monotone level 1 objective function if it can be written as

    f(Q)=g(Q, Q)

where  is some measure of central tendency and g is a function that is non-decreasing in both of its arguments. One can check that the impact measure  falls under this category of functions as do many of its variants.

The central observation for computing the tight optimistic estimator for monotone level 1 functions is that the optimum value must be attained on a sub-multiset that contains a consecutive segment of elements of Q from the top element w.r.t. y down to some cut-off element. Formally, let us define the top sequence of sub-multisets of Q as
T_i={y_m-i+1,…, y_m}
for i ∈m and note the following observation:

Let f be a monotone level 1 objective function.
Then the tight optimistic estimator of f can be computed as the maximum value on the top sequence, i.e.,
f(Q)=max{f(T_i)  i ∈m}.



Let R ⊆ Q be of size k with R={y_i_1,…,y_i_k}. Since y_i_j≤ y_m-j+1, we have for the top sequence element T_k that R  T_k and, hence, R≤T_k implying

    f(R)=g(k,R) ≤ g(k,T_l)=f(T_k)  .

It follows that for each sub-multiset of Q there is a top sequence element of at least equal objective value.

From this insight it is easy to derive an m algorithm for computing the tight optimistic estimator under the additional assumption
that we can compute g and the “incremental central tendency problem” (i,Q,(T_1,…,T_i-1) ↦T_i in constant time.
Note that computing the incremental problem in constant time implies to only access a constant number of target values and of the previously computed central tendency values. This can for instance be done for = via the incremental formula T_i=((i-1) T_i-1+y_m-i+1)/i or for = through direct index access of either of the two values y_m-⌊(i-1)/2 ⌋ or y_m-⌈ (i-1)/2 ⌉.
Since, according to Prop. <ref>, we have to evaluate f only for the m candidates T_i to find f(Q) we can do so in time m by solving the problem incrementally for i=1,…,m.
The same overall approach can be readily generalized for objective functions that are monotonically decreasing in the central tendency or those that can be written as the maximum of one monotonically increasing and one monotonically decreasing level 1 function. However, it breaks down for objective functions that depend on more than just size and central tendency—which inherently is the case when we want to incorporate dispersion-control.



 §.§ Dispersion-corrected objective functions based on the median


We will now extend the previous recipe for computing the tight optimistic estimator to objective functions that 
depend not only on subpopulation size and central tendency but also on the target value dispersion in the subgroup.
Specifically, we focus on the median as measure of central tendency and consider functions that are both monotonically increasing in the described subpopulation size and monotonically decreasing in some dispersion measure around the median.
To precisely describe this class of functions, we first have to formalize the notion of dispersion measure around the median.
For our purpose the following definition suffices.
Let us denote by Y the multiset of absolute differences to the median of a multiset Y ∈^, i.e., Y={y_1-Y,…,y_m-Y}.


We call a mapping ^ a dispersion measure around the median if Y is monotone with respect to the multiset of absolute differences to its median Y, i.e., if YZ then Y≤Z.

One can check that this definition contains the measures median absolute deviation around the median Y=Y, the root mean of squared deviations around the median Y={x^2  x ∈Y}^1/2, as well as the mean absolute deviation around the median Y=Y.[We work here with the given definition of dispersion measure because of its simplicity. Note, however, that all subsequent arguments can be extended in a straightforward way to a wider class of dispersion measures by considering the multisets of positive and negative deviations separately. This wider class also contains the interquartile range and certain asymmetric measures, which are not covered by Def. <ref>.]
Based on Def. <ref> we can specify the class of objective functions that we aim to tackle as follows: we call a function fP a dispersion-corrected or level 2 objective function (based on the median) if it can be written as

    f(Q)=g(Q,Q,Q)

where  is some dispersion measure around the median and g^3 is a real function that is non-decreasing in its first argument and non-increasing in its third argument (without any monotinicity requirement for the second argument).


Our recipe for optimizing these functions is then to consider only subpopulations R ⊆ Q that can be formed by selecting all individuals with a target value in some interval. Formally, for a fixed index z ∈{1,…,m} define m_z ≤ m as the maximal cardinality of a sub-multiset of the target values that has median index z, i.e.,

    m_z=min{2z,2(m-z)+1} .

Now, for k ∈m_z, let us define Q^k_z as the set with k consecutive elements around index z.
That is 

    Q^k_z={y_z-⌊k-1/2⌋, …, y_z, …, y_z+ ⌈k-1/2⌉} .

With this we can define the elements of the median sequence Q_z as those subsets of the form of Eq. (<ref>) that maximize f for some fixed index z ∈m.
That is, Q_z=Q^k^*_z_z where k^*_z ∈m_z is minimal with

    f(Q^k^*_z_z)=g(k^*_z, y_z,Q^k^*_z_z)=max{f(Q^k_z)  k ∈m_z} .

Thus, the number k^*_z is the smallest cardinality that maximizes the trade-off of size and dispersion encoded by g (given the fixed median y_z=Q^k_z for all k).
Fig. <ref> shows an exemplary median sequence based on 21 random target values.
In the following proposition we note that, as desired, searching the median sequence is sufficient for finding optimal subsets of Q.

Let f be a dispersion-corrected objective function based on the median.
Then the tight optimistic estimator of f can be computed as the maximum value on the median sequence, i.e.,
f(Q)=max{f(Q_z)  z ∈m}.



For a sub-multiset R ⊆ Q let us define the gap count γ(R) as 

    γ(R)={y ∈ Q ∖ R min R < y < max R} .

Let O ⊆ Q be an f-maximizer with minimal gap count, i.e., f(R) < f(O) for all R with γ(R)<γ(O). 
Assume that γ(O)>0. That means there is a y ∈ Q ∖ O such that min O < y < max O. Define

    S=
    
    (O ∖{min O}) ∪{y},    if y ≤O
    
    (O ∖{max O}) ∪{y},    otherwise .

Per definition we have S=O and S=O. Additionally, we can check that SO, and, hence, S≤Q. This implies that 

    f(S)=g(S,S,S)≥ g(O,O,O)=f(O)  .

However, per definition of S it also holds that γ(S) < γ(O), which contradicts that O is an f-optimizer with minimal gap count. Hence, any f-maximizer O must have a gap count of zero. In other words, O is of the form O=Q^k_z as in Eq. (<ref>) for some median z ∈m and some cardinality k ∈m_z and per definition we have f(Q_z) ≥ f(O) as required. 
Consequently, we can compute the tight optimistic estimator for any dispersion-corrected objective function based on the median in time m^2 for subpopulations of size m—again, given a suitable incremental formula for .
While this is not generally a practical algorithm in itself, it is a useful departure point for designing one.
In the next section we show how it can be brought down to linear time when we introduce some additional constraints on the objective function.




 §.§ Reaching linear time—objectives based on dispersion-corrected coverage


Equipped with the general concept of the median sequence, we can now address the special case of 
dispersion-corrected objective functions where the trade-off between the subpopulation size and target value dispersion is captured by a linear function of size and the sum of absolute differences from the median. Concretely, let us define the dispersion-corrected coverage (w.r.t. absolute median deviation) by

    Q=Q/P(1-Q/P)_+ =  ( Q/P-Q/P)_+

where Q=∑_y ∈ Qy-Q denotes the sum of absolute deviations from the median. We then consider objective functions based on the dispersion-corrected coverage of the form

    f(Q)=g(Q,Q)

where g is non-decreasing in its first argument. Let us note, however, that we could replace the  function by any linear function that depends positively on Q and negatively on . It is easy to verify that function of this form also obey the more general definition of level-2 objective functions given in Sec. <ref>, and, hence can be optimized via the median sequence.

The key to computing the tight optimistic estimator f in linear time for functions based on dispersion-corrected coverage is then that the members of the median sequence Q_z can be computed incrementally in constant time. Indeed, we can prove the following theorem, which states that the optimal size for a multiset around median index z is within 3 of the optimal size for a multiset around median index z+1—a fact that can also be observed in the example given in Fig. <ref>.


Let f be of the form of Eq. (<ref>). For z ∈m-1 it holds for the size k_z^* of the f-optimal multiset with median z that 

    k_z^* ∈{max(0,k^*_z+1-3),…,min(m_z,k^*_z+1+3)} .


One idea to prove this theorem is to show that a) the gain in f for increasing the multiset around a median index z is alternating between two discrete concave functions and b) that the gains for growing multisets between two consecutive median indices are bounding each other.
For an intuitive understanding of this argument, Fig. <ref> shows for four different median indices z ∈{10,11,12,13} the dispersion-corrected coverage for the sets Q^k_z as a function in k. On closer inspection, we can observe that when considering only every second segment of each function graph, the corresponding -values have a concave shape.
A detailed proof, which is rather long and partially technical, can be found in Appendix <ref>.


It follows that, after computing the objective value of Q_m trivially as f(Q_m)=g(1/P,y_m), we can obtain f(Q_z-1) for z=m,…,2 by checking the at most seven candidate set sizes given by Eq. (<ref>) as

    f(Q^z-1)=max{f(Q^k^-_z_z-1),…, f(Q^k_z^+_z-1)}

with k_z^-=max(k_z^*-3,1) and k_z^+=min(k_z^*+3,m_z).
It remains to see that we can compute individual evaluations of f in constant time (after some initial m pre-processing step). As a general data structure for quickly computing sums of absolute deviations from a center point, we can define for i ∈m the left error e_l(i) and the right error e_r(i) as

    e_l(i)=∑_j=1^i-1 y_i-y_j,  e_r(i)=∑_j=i+1^m y_j-y_i  .

Note that we can compute these error terms for all i ∈m in time m via the recursions

    e_l(i)   =e_l(i-1)+(i-1)(y_i-y_i-1)
    
    e_r(i)   =e_r(i+1)+(m-i)(y_i+1-y_i)

and e_l(1)=e_r(m)=0. Subsequently, we can compute sums of deviations from center points of arbitrary subpopulations in constant time, as the following statement shows (see Appendix <ref> for a proof).


Let Q={y_1,…,y_a, …, y_z, …, y_b,…, y_m} be a multiset with 1 ≤ a<z<b ≤ m and y_i≤ y_j for i≤ j. Then the sum of absolute deviations to y_i of all elements of the submultiset {y_a, …, y_z, …, y_b} can be expressed as

    ∑_i=a^b y_z-y_i   =e_l(z)-e_l(a)-(a-1)(y_z-y_a)
       +e_r(z)-e_r(b)-(m-b)(y_b-y_z)  .


With this we can compute k ↦ f(Q_z^k) in constant time (assuming g can be computed in constant time). Together with Prop. <ref> and Thm. <ref> this results in a linear time algorithm for computing Q ↦f(Q) (see Alg. <ref> for a pseudo-code that summarizes all ideas).




§ DISPERSION-CORRECTED SUBGROUP DISCOVERY IN PRACTICE









[ph!]


!
7cDataset     6cSelection Bias     5cEfficiency

(lr)1-7 (lr)8-13 (lr)14-18
       Name            Target     P     Π     P     P                                          a_eff               t_0                   t_1               
 1       abalone         rings     4,177          69               9                2.359          0.544     0.191                11                      11                      2.257                      1.662        1                  848,258            690,177            304          339               

2       ailerons        goal     13,750         357              -0.0008        0.000303     0.906     0.59                 -0.0007               -0.0006      0.000288                 0.000198     0.3                1,069,456          54,103             6,542                 460      

3       autoMPG8        mpg     392            24               22.5             6.524          0.497                0.497                29                      29                      4.791                      4.791                   1                  96                   67                   0.11                  0.09     

4       baseball        salary     337            24               740              954.386        0.362     0.003                1550                    2500         1245.092     0            1                  117                  117                  0.22                  0.21     

5       california      med. h. value     20,640       72               179,700        88,354       0.385     0.019              262,500               500,001      94261        294,00       0.4                1,368,662          65,707             2,676                 368      

6       compactiv       usr     8,192          202              89               9.661          0.464                0.603     94           93                      7.8                        3.472        0.5                2,458,105          59,053             5,161                 208      

7       concrete        compr. strength     1,030          70               34.4             13.427         0.284     0.1291               48.97                   50.7         12.744                     9.512        1                  512,195            221,322            43.9                  35.8     

8       dee             consume     365            60               2.787            0.831          0.523     0.381              3.815                 4.008        0.721                    0.434        1                  18,663             2,653              2.05                  1.29     

9       delta_ail      sa     7,129        66               -0.0001        0.000231     0.902     0.392                0.0001                0.0002       0.000226                 0.000119     1                  45,194             2,632              33.3                  6.11     

10      delta_elv      se     9,517        66               0.001            0.00198        0.384     0.369                0.002                   0.002                 0.00112                    0.00108      1                  10145                1,415              8.9                   4.01     

11      elevators       goal     16,599       155              0.02             0.00411        0.113                0.283     0.03         0.021                 0.00813      0.00373      0.05               6,356,465          526,114            13,712                2,891    

12      forestfires     area     517            70               0.52             12.832         0.01      0.002              86.45                 278.53       56.027       0            1                  340,426            264,207            23           23.7              

13      friedman        output     1,200        48               14.651           4.234          0.387     0.294                18.934                  19.727       3.065                      2.73         1                  19,209             2,489              3.23                  1.56     

14      house           price     22,784       160              33,200         28,456         0.56                 0.723     45,200       34,000                40,576       27,214       0.002              1,221,696          114,566            7,937                 1,308    

15      laser           output     993          42               46               35.561         0.32      0.093                109                     135          40.313       15.662       1                  2,008              815                0.96                  0.83     

16      mortgage        30 y. rate     1,049        128              6.71             2.373          0.256     0.097                11.61                   14.41        2.081                      0.98         1                  40,753             1,270              11.6                  1.59     

17      mv              y     40,768       79               -5.02086         8.509          0.497     0.349              0.076                   0.193        8.541        2.032        1                  6,513              1,017              31.9                  13.2     

18      pole            output     14,998       260              0                28.949         0.40      0.24                 100                     100                     38.995       16.692       0.2                1,041,146          2,966              2,638                 15       

19      puma32h         thetadd6     8,192        318              0.000261       0.023          0.299     0.244                0.026                   0.031        0.018                      0.017        0.4                3,141,046          5,782              2,648                 15.5     

20      stock           company10     950            80               46.625           5.47           0.471     0.337                52.5                    54.375       3.741                      2.515        1                  85,692             1,822              12.5                  1.56     

21      treasury        1 m. def. rate     1,049        128              6.61             2.473          0.182                0.339     13.16        8.65                  2.591        0.863        1                  49,197             9,247              14.8                  5.91     

22      wankara         mean temp.     321        87               47.7             12.753         0.545     0.296                60.6                    67.6         8.873                      4.752        1                  191,053            4,081              11.9                  1.24     

23      wizmir          mean temp.     1,461        82               60               12.622         0.6       0.349                72.9                    78.5         8.527                      3.889        1                  177,768            1,409              38.5                  1.48     
 24      binaries        delta E     82             499              0.106            0.277          0.305              0.378     0.43         0.202                   0.373        0.118        0.5                4,712,128          204                1,200                 0.29     

25      gold            Evdw-Evdw0     12,200       250              0.131            0.088          0.765     0.34                 0.217                   0.234        0.081                      0.0278       0.4                1,498,185          451                5,650                 3.96     



Datasets with corresponding population size (P), number of base propositions (Π), global median (P) and mean absolute median deviation (P) followed by coverage (, ), median (, ), and mean absolute median deviation (, ) for best subgroup w.r.t. non-dispersion corrected function f_0 and dispersion-corrected function f_1, respectively; bold-face indicates higher coverage and median and lower dispersion, underlines indicate higher dispersion than in global population; final column segment contains accuracy parameter used in the efficiency study (a_eff) as well as number of expanded nodes (, ) and computation time in seconds (t_0, t_1) for optimistic estimator based on top sequence f_0 and tight optimistic estimator f_1, respectively—in both cases when optimizing f_1; depth-limit of 10 is used for all datasets with a<1, no depth-limit otherwise.


The overall result of Sec. <ref> is an efficient algorithm for dispersion-corrected subgroup discovery which, e.g., allows us to replace the coverage term in standard objective functions by the dispersion-corrected coverage.
To evaluate this efficiency claim as well as the value of dispersion-correction, let us consider as objective the normalized and dispersion-corrected impact function based on the median, i.e., f_1(Q)=QQ where  is the positive relative median shift 

    Q= ( Q-P/P-P )_+  .

This function obeys Eq. (<ref>); thus, its tight optimistic estimator can be computed using the linear time algorithm from Sec. <ref>.

The following empirical results were gathered by applying it to a range of publicly available real-world datasets.[Datasets contain all regression datasets from the KEEL repository <cit.> with at least 5 attributes and two materials datasets from the Nomad Repository <nomad-coe.eu/>; see Tab. <ref>. Implementation available in open source Java library realKD <bitbucket.org/realKD/>. Computation times determined on MacBook Pro 3.1 GHz Intel Core i7.]
We will first investigate the effect of dispersion-correction on the output before turning to the effect of the tight optimistic estimator on the computation time. 



 §.§ Selection Bias of Dispersion-Correction and its Statistical Merit


To investigate the selection bias of f_1 let us also consider the non-dispersion corrected variant f_0(Q)=QQ where we simply replace the dispersion-corrected coverage by the ordinary coverage. This function is a monotone level 1 function, hence, its tight optimistic estimator f_0 can be computed in linear time using the top sequence approach.
Fig. <ref> shows the characteristics of the optimal subgroups that are discovered with respect to both of these objective functions (see also Tab. <ref> for exact values) where for all datasets the language of closed conjunctions  has been used as description language.

The first observation is that—as enforced by design—for all datasets the mean absolute deviation from the median is lower for the dispersion-corrected variant (except in one case where both functions yield the same subgroup). On average the dispersion for f_1 is 49 percent of the global dispersion, whereas it is 113 percent for f_0, i.e., when not optimizing the dispersion it is on average higher in the subgroups than in the global population.
When it comes to the other subgroup characteristics, coverage and median target value, the global picture is that f_1 discovers somewhat more specific groups (mean coverage 0.3 versus 0.44 for f_0) with higher median shift (on average 0.73 normalized median deviations higher). 
However, in contrast to dispersion, the behavior for median shift and coverage varies across the datasets.
In Fig. <ref>, the datasets are ordered according to the difference in subgroup medians between the optimal subgroups w.r.t. f_0 and those w.r.t. f_1.










This ordering reveals the following categorization of outcomes: When our description language is not able to reduce the error of subgroups with very high median value, f_1 settles for more coherent groups with a less extreme but still outstanding central tendency.
On the other end of the scale, when no coherent groups with moderate size and median shift can be identified, the dispersion-corrected objective selects very small groups with the most extreme target values.
The majority of datasets obey the global trend of dispersion-correction leading to somewhat more specific subgroups with higher median that are, as intended, more coherent.






To determine based on these empirical observations, whether we should generally favor dispersion correction, we have to specify an application context that specifies the relative importance of coverage, central tendency, and dispersion. For that let us consider the common statistical setting in which we do not observe the full global population P but instead subgroup discovery is performed only on an i.i.d. sample P' ⊆ P yielding subpopulations Q'=σ(P'). While σ has been optimized w.r.t. the statistics on that sample Q' we are actually interested in the properties of the full subpopulation Q = σ(P).
For instance, a natural question is what is the minimal y-value that we expect to see in a random individual q ∈ Q with high confidence.
That is, we prefer subgroups with an as high as possible threshold  such that a random q ∈ Q satisfies with probability[The probability is w.r.t. to the  distribution with which the sample P' ⊆ P is drawn.] 1-δ that y(q) ≥.
This criterion gives rise to a natural trade-off between the three evaluation metrics through the empirical Chebycheff inequality <cit.>, according to which we can compute such a value as Q'-ϵ(Q') where

    ϵ(Q') = √((Q'^2-1)(Q')/Q'^2δ - Q')

and Y=∑_y ∈ Y (y-Y)^2 / (Y-1) is the sample variance. Note that this expression is only defined for sample subpopulations with a size of at least 1/δ. For smaller subgroups our best guess for a threshold value would be the one derived from the global sample P'-ϵ(P') (which we assume to be large enough to determine an ϵ-value). This gives rise to the following standardized lower confidence bound score  that evaluates how much a subgroup improves over the global  value:

    Q'=( Q'-P'/√((P')) )_+ whereQ'=Q'-ϵ(Q')    , if ϵ(Q')  defined
    P'-ϵ(P')    , otherwise .


The plot on the left side of Fig. <ref> shows the score values of the optimal subgroup w.r.t. to f_1 (_1) and f_0 (_0) using confidence parameter δ=0.05. Except for three exceptions (datasets 3,4, and 12), the subgroup resulting from f_1 provides a higher lower bound than those from the non-dispersion corrected variant f_0. That is, the data shows a strong advantage for dispersion correction when we are interested in selectors that mostly select individuals with a high target value from the underlying population P.
In order to test the significance of these results, we can employ the Bayesian sign-test <cit.>, a modern alternative to classic frequentist null hypothesis tests that avoids many of the well-known disadvantages of those <cit.>. With Bayesian hypothesis tests, we can directly evaluate the posterior probabilities of hypotheses given our experimental data instead of just rejecting a null hypothesis based on some arbitrary significance level. Moreover, we differentiate between sample size and effect size by the introduction of a region of practical equivalence (rope). Here, we are interested in the relative difference z̃=(_1 - _0)/(max{_0, _1}) on average for random subgroup discovery problems.

Using a conservative choice for the rope, we call the two objective functions practically equivalent if the mean z̃-value is at most r=0.1. Choosing the prior belief that f_0 is superior, i.e., z̃ < -r, with a prior weight of 1, the procedure yields based on our 25 test datasets the posterior probability of approximately 1 that z̃ > r on average (see the right part of Fig. <ref> for in illustration of the posterior belief). Hence, we can conclude that dispersion-correction improves the relative lower confidence bound of target values on average by more than 10 percent when compared to the non-dispersion-corrected function.




 §.§ Efficiency of the Tight Optimistic Estimator

To study the effect of the tight optimistic estimator, let us compare its performance to that of a baseline estimator that can be computed with the standard top sequence approach. Since f_1 is upper bounded by f_0, f_0 is a valid, albeit non-tight, optimistic estimator for f_1 and can thus be used for this purpose.


The exact speed-up factor is determined by the ratio of enumerated nodes for both variants as well as the ratio of computation times for an individual optimistic estimator computation.
While both factors determine the practically relevant outcome, the number of nodes evaluated is a much more stable quantity, which indicates the full underlying speed-up potential independent of implementation details. 
Similarly, “number of nodes evaluated” is also an insightful unit of time for measuring optimization progress.
Therefore, in addition to the computation time in seconds t_0 and t_1, let us denote by , ⊆ the set of nodes enumerated by branch-and-bound using f_0 and f_1, respectively—but in both cases for optimizing the dispersion-corrected objective f_1.
Moreover, when running branch-and-bound with optimistic estimator f_i, let us denote by σ^*_i(n) and σ^+_i(n) the best selector found and the top element of the priority queue (w.r.t. f_i), respectively, after n nodes have been enumerated.

The plot on the left side of Fig. <ref> shows the speed-up factor t_1/t_0 on a logarithmic axis for all datasets in increasing order along with the potential speed-up factors / (see Tab. <ref> for numerical values). There are seven datasets for which the speed-up turns out to be minor followed by four datasets with a modest speed-up factor of 2. For the remaining 14 datasets, however, we have solid speed-up factors between 4 and 20 and in four cases immense values between 100 and 4,000.
This demonstrates the decisive potential effect of tight value estimation even when compared to another non-trivial estimator like f_0 (which itself improves over simpler options by orders of magnitude; see ). 
Similar to the results in Sec. <ref>, the Bayesian sign-test for the normalized difference z=(t_1-t_0)/max{t_1,t_0} with the prior set to practical equivalence (z ∈ [-0.1,0.1]) reveals that the posterior probability of f_1 being superior to f_0 is approximately 1.


In almost all cases the potential speed-up given by the ratio of enumerated nodes is considerably higher than the actual speed-up, which shows
that, despite the same asymptotic time complexity, an individual computation of the tight optimistic estimator is slower than the simpler top sequence based estimator—but also indicates that there is room for improvements in the implementation.
When zooming in on the optimization progress over time for the binaries dataset, which exhibits the most extreme speed-up (right plot in Fig. <ref>),
we can see that not only does the tight optimistic estimator close the gap between best current selector and current highest potential selector much faster—thus creating the huge speed-up factor—but also that it causes better solutions to be found earlier.
This is an important property when we want to use the algorithm as an anytime algorithm, i.e., when allowing the user to terminate computation preemptively, which is important in interactive data analysis systems.
This is an advantage enabled specifically by using the tight optimistic estimators in conjunction with the best-first node expansion strategy.



§ CONCLUSION

During the preceding sections, we developed and evaluated an effective algorithm for simultaneously optimizing size, central tendency, and dispersion in subgroup discovery with a numerical target.
This algorithm is based on two central results: 1) the tight optimistic estimator for any objective function that is based on some dispersion measure around the median can be computed as the function's maximum on a linear-sized sequence of sets—the median sequence (Prop. <ref>); and 2) for objective functions based on the concept of the dispersion-corrected coverage w.r.t. the absolute deviation from the median, the individual sets of the median sequence can be generated in incremental constant time (Thm. <ref>).



Among the possible applications of the proposed approach, the perhaps most important one is to replace the standard coverage term in classic objective functions by the dispersion-corrected coverage, i.e., the relative subgroup size minus the relative subgroup dispersion, to reduce the error of result subgroups—where error refers to the descriptive or predictive inaccuracy incurred when assuming the median value of a subgroup for all its members.
As we saw empirically for the impact function (based on the median), this correction also has a statistical advantage resulting in subgroups where we can assume greater target values for unseen group members with high confidence.
In addition to enabling dispersion-correction to known objective functions, the presented algorithm also provides novel degrees of freedom, which might be interesting to exploit in their own right:
The dependence on the median is not required to be monotone, which allows to incorporate a more sophisticated influence of the central tendency value than simple monotone average shifts. For instance, given a suitable statistical model for the global distribution, the effect of the median could be a function of the probability ℙ[Q], e.g., its Shannon information content.
Furthermore, the feasible dispersion measures allow for interesting weighting schemes, which include possibilities of asymmetric effects of the error (e.g., for only punishing one-sided deviation from the median).





Regarding the limitations of the presented approach, let us note that it cannot be directly applied to the previously proposed dispersion-aware functions, i.e., the t-score 𝚝𝚜𝚌(Q)=√(Q)(Q-P)/Q and the mmad score for ranked data 𝚖𝚖𝚍(Q)=Q/(2Q+Q). 
While both of these functions can be optimized via the median sequence approach (assuming a t-score variant based on the median), we are lacking an efficient incremental formula for computing the individual function values for all median sequence sets, i.e., a replacement for Thm. <ref>.
Though finding such a replacement in future research is conceivable, 
this leaves us for the moment with a quadratic time algorithm (in the subgroup size) for the tight optimistic estimator, which is not generally feasible (although potentially useful for smaller datasets or as part of a hybrid optimistic estimator, which uses the approach for sufficiently small subgroups only).


Since they share basic monotonicities, it is possible to use functions based on dispersion-corrected coverage as an optimization proxy for the above mentioned objectives. For instance, the ranking of the top 20 subgroups w.r.t. the dispersion-corrected binomial quality function, 𝚍𝚌𝚋(Q)=√(Q)(Q-P), turns out to have a mean Spearman rank correlation coefficient with the median-based t-score of apx. 0.783 on five randomly selected test datasets (delta_elv, laser, stock, treasury, gold). However, a more systematic understanding of the differences and commonalities of these functions is necessary to reliably replace them with one another. Moreover, the correlation deteriorates quite sharply when we compare to the original mean/variance based t-score (mean Spearman correlation coefficient 0.567), which points to the perhaps more fundamental limitation of the presented approach for dispersion-correction: it relies on using the median as measure of central tendency.
While the median and the mean absolute deviation from the median are an interpretable, robust, and sound combination of measures (the median of a set of values minimizes the sum of absolute deviations), the mean and the variance are just as sound, are potentially more relevant when sensitivity to outliers is required, and provide a wealth of statistical tools (e.g., the empirical Chebyshev's inequality used above).


Hence, a straightforward but valuable direction for future work is the extension of efficient tight optimistic estimator computation to dispersion-correction based on the mean and variance. A basic observation for this task is that objective functions based on dispersion measures around the mean must also attain their maximum on gap-free intervals of target values. However, for a given collection of target values, there is a quadratic number of intervals such that a further idea is required in order to attain an efficient, i.e., (log-)linear time algorithm.
Another valuable direction for future research is the extension of consistency and error optimization to the case of multidimensional target variables where subgroup parameters can represent complex statistical models <cit.>. While this setting is algorithmically more challenging than the univariate case covered here, the underlying motivation remains: balancing group size and exceptionality, i.e., distance of local to global model parameters, with consistency, i.e., local model fit, should lead to the discovery of more meaningful statements about the data and the underlying domain.




  
Acknowledgements
The authors thank the anonymous reviewers for their useful and constructive suggestions.
Jilles Vreeken and Mario Boley are supported by the Cluster of Excellence “Multimodal Computing and Interaction” within the Excellence Initiative of the German Federal Government. 
Bryan R. Goldsmith acknowledges support from the Alexander von Humboldt-Foundation with a Postdoctoral Fellowship.

Additionally, this work was supported through the European Union's Horizon 2020 research and innovation program under
grant agreement no. 676580 with The Novel Materials Discovery (NOMAD) Laboratory, a European Center of
Excellence.

spbasic 



§ PROOF OF THEOREM <REF>


In order to proof Thm. <ref>, let us start by noting that for functions of the form of Eq. (<ref>), finding the set size k^*_z corresponds to maximizing the dispersion-corrected coverage among all multisets with consecutive elements around median y_z (as defined in Eq. <ref>).
In order to analyze this problem, let us write

    h_z(k)=Q^k_z=Q^k_z/P - Q^k_z/P

for the dispersion-corrected coverage of the multiset Q^k_z.
Let Δ h_zm_z denote the difference or gain function of h_z, i.e., Δ h_z(k)=h_z(k)-h_z(k-1) where we consider Q^0_z=∅ and, hence, h_z(0)=0. With this definition we can show that h_z is alternating between two concave functions, i.e., considering either only the even or only the odd subset of its domain, the gains are monotonically decreasing. More precisely:


For all k ∈m_z∖{1,2} we have that Δ h_z(k) ≤Δ h_z(k-2).


For k ∈m_z, let us denote by q^k_z the additional y-value that Q_z^k contains compared to Q_z^k-1 (considering Q_z^0=∅), i.e., Q_z^k ∖ Q_z^k-1={q_z^k}. We can check that 

    q_z^k=
    
    q_z-⌊k-1/2⌋,    k  odd
    
    q_z+⌈k-1/2⌉,    k  even .

With this and using the shorthands n=P and d=P we can write

    Δ h_z(k) - Δ h_z(k-2)
    = h_z(k)-h_z(k-1)-(h_z(k-2)-h_z(k-3))
    
    =   k/n-Q^k_z/d-k-1/n+Q^k-1_z/d-k-2/n+Q^k-2_z/d+k-3/n-Q^k-3_z/d
    
    =   1/n (k-k+1-k+2+k-3  )_=0+1/d (Q_z^k-2- Q^k_z+Q_z^k-1-Q_z^k-3 )
    
    =   1/d( -q_z^k-y_z - q^k-1_z-y_z + q_z^k-1-y_z + q_z^k-2-y_z)
    
    =   1/d( -q_z^k-y_z + q_z^k-2-y_z)
    case k odd
    =   1/d( -(y_z-y_z-⌊k-1/2⌋) + (y_z-y_z-⌊k-3/2⌋) ) = y_z-⌊k-1/2⌋ - y_z-⌊k-3/2⌋≤ 0
    case k even
    =   1/d( -(y_z+⌈k-1/2⌉-y_z) + (y_z+⌈k-3/2⌉-y_z) ) = y_z+⌈k-3/2⌉ - y_z+⌈k-1/2⌉≤ 0


One important consequence of this fact is that the operation of growing a set around median z by two elements—one to the left and one to the right—has monotonically decreasing gains.
In other words, the smoothed function h_z(k)=h_z(k)+h_z(k-1) is concave or formally

    Δ h_z(k) + Δ h_z(k-1) ≥Δ h_z(k+1)+Δ h_z(k)  .

Moreover, we can relate the gain functions of consecutive median indices as follows.


Let z ∈m∖{1} and k ∈m_z-1∖{1,2,3}. It holds that

    Δ h_z-1(k-2)+Δ h_z-1(k-3)    ≥Δ h_z(k)+Δ h_z(k-1)
    Δ h_z-1(k)+Δ h_z-1(k-1)    ≤Δ h_z(k-2)+Δ h_z(k-3)







For this proof, let us use the same shorthands as in the proof of Lemma <ref> and start by noting that for all i ∈m and k ∈m_z∖{1} we have the equality

    Δ h_i(k)+Δ h_i(k-1) = 2/n-q_i^k-y_i+q_i^k-1-y_i/d

which we can see by extending

    Δ h_i(k)+Δ h_i(k-1)    =h_i(k)-h_i(k-1)+h_i(k-1)-h_i(k-2) 
       =k-k+2/n-Q^k_i-Q^k-2_i/d=2/n-q^k_i-y_i+q^k-1_i-y_i/d .

We can then show Eq. (<ref>) by applying Eq. (<ref>) two times to

    Δ h_z-1(k-2)+Δ h_z-1(k-3) - (Δ h_z(k)+Δ h_z(k-1))
    
    
    =   1/d ( -q^k-2_z-1-y_z-1 - q^k-3_z-1-y_z-1 + q^k_z-y_z + q^k-1_z-y_z )
    and finally by checking separately the case k odd
    =   1/d ( y_z-1-k-3/2-y_z-1+y_z-1 - y_z-1+k-4/2 + y_z - y_z-k-1/2 + y_z+k-2/2⌉-y_z )
    
    =   1/d ( y_z-k-1/2 - y_z-k-1/2_=0  + y_z-1+k/2 - y_z-3+k/2_≥ 0 ) ≥ 0
    and the case k even
    =   1/d(y_z-1- y_z-1+k-3/2 +  y_z-1-k-4/2-y_z-1 + y_z+k-1/2 - y_z + y_z - y_z - k-2/2)
    
    =   1/d(y_z+1-k/2 - y_z +1 - k/2_=0 + y_z+k-1/2 - y_z-2+k-1/2_≥ 0) ≥ 0  .

Similarly, for Eq. (<ref>) by applying Eq. (<ref>) two times we can write

    Δ h_z-1(k)+Δ h_z-1(k-1) - (Δ h_z(k-2) + Δ h_z(k-3))
    
    =   1/d( -q_z-1^k-y_z-1 - q_z-1^k-1 - y_z-1 + q_z^k-2-y_z + q_z^k-3-y_z)
    
    
    =   1/d( y_z-1-k-1/2 - y_z+1-k-1/2_≤ 0 + y_z-2+k/2-y_z-2+k/2_=0) ≤ 0,    k  odd
    1/d( y_z-k/2 - y_z+2-k/2_≤ 0 + y_z-1+k-1/2 - y_z-1+k-1/2_=0) ≤ 0,    k  even


Combining all of the above we can finally proof our main result as follows.

We start by showing that every k ∈m_z+1 with k<k^*_z-3 can not be an optimizer of h_z+1. It follows that k^*_z-3 ≤ k^*_z+1, and, hence, k^*_z≤ k^*_z+1+3 as required for the upper bound. Indeed, we have

    h_z+1(k)    = h_z+1(k+2) - (Δ h_z+1(k+2)+Δ h_z+1(k+1))
       ≤ h_z+1(k+2) - (Δ h_z+1(k^*_z-2) + Δ h_z+1(k^*_z-3))    (by Eq. (<ref>))
       ≤ h_z+1(k+2) -  (Δ h_z(k^*_z) + Δ h_z(k^*_z-1))_>0  by def. of k^*_z  <h_z+1(k+2)  .    (by Lm. <ref>)

Analogously, for the lower bound, we show that every k ∈m_z+1 with k>k^*_z+3 can not be the smallest optimizer of h_z+1. It follows that k^*_z+3 ≥ k^*_z+1, and, hence, k^*_z≥ k^*_z+1-3 as required. Indeed, we can write

    h_z+1(k)    = h_z+1(k-2)+Δ h_z+1(k)+Δ h_k+1(k-1) 
       ≤ h_z+1(k-2) + Δ h_z+1(k^*_z+4) + Δ h_z+1(k^*_z+3)     (by Eq. (<ref>))
       ≤ h_z+1(k-2) + Δ h_z(k^*_z+2) + Δ h_z(k^*_z+1)_≤ 0  by def. of k^*_z≤ h_z+1(k-2)    (by Lm. <ref>)





§ ADDITIONAL PROOFS



Using d_ij as a shorthand for y_j-y_i for i,j ∈m with i≤ j we can write

    e_l(z)-e_l(a)-(a-1)(y_z-y_a)+e_r(z)-e_r(b)-(m-b)d_zb
    
    =   ∑^z-1_i=1d_iz-∑^a-1_i=1d_ia-(a-1)d_az+∑_i=z+1^m d_zi-∑_i=b+1^md_bi-(m-b)d_zb
    
    =   ∑^z-1_i=ad_iz+∑^a-1_i=1(d_iz-d_ia)_d_az-(a-1)d_az+∑^b_i=z+1d_zi+∑^m_i=b+1(d_zi-d_bi)_d_zb-(m-b)d_zb
    
    =   ∑^z-1_i=ad_iz+(a-1)d_az-(a-1)d_az+∑^b_i=z+1d_zi+(m-b)d_zb-(m-b)d_zb
    
    =   ∑^z-1_i=ad_iz + ∑^b_i=z+1d_zi=∑_i=a^b y_z-y_i





§ SUMMARY OF USED NOTATIONS




Symbol     Meaning     Defined in

·     cardinality of a set or absolute value of a number     - 

k     set of integers {1,…,k}     <ref>

(x)_+     max{x,0} for a real-valued expression x     <ref>

     element-wise less-or-equal relation for multisets of real values     <ref>

X     power set of a set X, i.e., set of all of its subsets     - 

^X     set of all multisets containing elements from set X     - 

σ, φ     subgroup selectors σ, φP{,}     <ref>

     measure of central tendency     <ref>

     measure of dispersion     <ref>

e_l(i), e_r(i)     left and right cumulative errors of target values up to value i     <ref>

f     objective function     <ref>

f     tight optimistic estimator of objective function f     <ref>

m     number of elements in subpopulation Q     <ref>

m_z     maximal size parameter k for consecutive value set Q_z^k     <ref>

k^*_z     f-maximizing size parameter k for consecutive value set Q_z^k     <ref>

y     numeric target attribute yP     <ref>

y_i     i-th target value of subpopulation w.r.t. ascending order     <ref>

P     global population of given subgroup discovery problem     <ref>

Q     some subpopulation Q ⊆ P     <ref>

Q_z     median sequence element with median index z     <ref>

Q^k_z     submultiset of Q with k consecutive elements around index z     <ref>

T_i     top sequence element i, i.e., T_i={y_m-i+1,…, y_m}     <ref>

Y     real-valued multiset     <ref>

Y     multiset of differences of elements in Y to its median     <ref>

,      description language and language of conjunctions     <ref>

     language of closed conjunctions     <ref>

Q     mean absolute deviation of y-values in Q to their median     <ref>

Q     coverage, i.e., relative size Q/P of subpopulation Q     <ref>

Q     dispersion-corrected coverage of subpopulation Q     <ref>

Q     arithmetic mean of y-values in Q     -

Q     impact, i.e., weighted mean-shift, of subpopulation Q     <ref>

Q     median of y-values in Q     <ref>

Q     sum of absolute deviations of y-values in Q to their median     <ref>





