PED-ANOVA: Efficiently Quantifying Hyperparameter Importance in Arbitrary Subspaces

1 Introduction

Following on the heels of widespread adoption of deep learning models in various industries and areas of research, Hyperparameter Optimization (HPO) [Bergstra and Bengio, 2012; Snoek et al., 2012; Bergstra et al., 2011; Lindauer et al., 2022; Watanabe, 2023] for deep learning has gained increasing prominence as the path forward for making deep learning more accessible and robust. In particular, recent research has highlighted the role that good hyperparameter (HP) space design can play in training strong models [Chen et al., 2018; Melis et al., 2018; Henderson et al., 2018]. In practice, while a large search space is necessary to find high-performance models [Zimmer et al., 2021], a reduced search space that retains only important HPs is essential for efficiently finding them [Perrone et al., 2019]. Therefore, it is crucial to understand the role that different HPs play in a search space.

This is the driving force behind previous research into the quantification of HP Importance (HPI) [Hutter et al., 2014;

Biedenkapp et al., 2017], which still remains a largely understudied section of HPO research. Several HPO frameworks [Biedenkapp et al., 2018; Akiba et al., 2019; Sass et al., 2022] have previously utilized functional ANOVA (fANOVA) [Hooker, 2007; Hutter et al., 2014] to provide a better interpretation of the role of different HPs, but the original f-ANOVA formulation is not very practical for the interpretation of specific subspaces of a search space. Such subspaces are often of particular interest to algorithm developers due to various properties, for example, the "local space" visualized in Figure 1 could represent a region of high performance. Nevertheless, prior works [Hutter et al., 2014; Biedenkapp et al., 2018] have attempted to overcome this and quantify HPI in specific subspaces using f-ANOVA. However, since their formulation did not constrain the calculations to subspaces of high interest, we argue that the results are biased towards unimportant subspaces. At the same time, obtaining an unbiased quantification of HPI in specific subspaces of interest, which we refer to as local spaces in contrast to the full global space, is mathematically non-trivial.

To overcome this issue, we first formally define local HPI as HPI in a local space and we derive a novel formulation of f-ANOVA to compute local HPI for arbitrary local spaces. Still, our formulation would require Monte-Carlo sampling

(b) An empirical γ-set PDF using KDE of f(x) = sin x

• (a) The true γ-set of f(x1, x2) = x21 + x22

Figure 2: Conceptual visualizations of the γ-set and the γ-set PDF. Left: the true γ-set in a 2D example. Darker is better in the figure and we colored the top-10% in red, which is the γ-set X γ (γ = 0.1) in this example. Right: the true γ-set and the γ-set PDF in a 1D example. The green dotted line shows the γ-quantile value f γ , which achieves the top-20% in this example, and the red solid lines are the γ-set X γ (γ = 0.2). The red dotted spaces below show the true γ-set PDF p(x|X γ ), but since we do not have the analytical form in practice, this PDF is estimated by KDE (the blue solid line).

in general and it is computationally intractable. Therefore, we show that local HPI is tractable without a Monte-Carlo sampling under some constraints and propose an algorithm that uses Pearson divergence (PED, [Pearson, 1900]) to enable a closed-form computation of HPI. In a series of experiments, we first verify that our algorithm correctly provides global and local HPI in a toy function. Then we demonstrate that our algorithm takes only less than a second for 105 data points while the prior f-ANOVA [Hutter et al., 2014] would take more than a week.

To provide a solid picture of how to use our method, we perform analysis on JAHS-Bench-201 [Bansal et al., 2022], which has one of the largest search spaces among HPO benchmarks. In the analysis, we find that it is suboptimal to design a search space relying only on global HPI because we potentially miss important HPs in a local space if the global HPI of these HPs are dominated by the most important HP. We demonstrate that local HPI plays a crucial role to avoid this issue. Furthermore, our method has several other possible applications such as (1) post-hoc analysis for HPO, (2) adaptive (e.g., meta-learned) search space reductions for faster HPO, and (3) exploratory data analysis. We discuss these in more detail in Appendix E, along with the advantages and limitations of our method.

In summary, the contributions of this paper are to:

1. reformulate local HPI mathematically and derive the general formula of local HPI,
2. provide a closed-form calculation for local HPI using PED that handles even 108 data points in a minute, and
3. benchmark performance compared with the original fANOVA.

To facilitate reproducibility, our implementation is available at https://github.com/nabenabe0928/local-anova/.

2 Background &RelatedWork

2.1 Preliminaries

Throughout this paper, we use the following terms:

1. γ-quantile value f γ : The function value f γ ∈ R that achieves the top-γ quantile with respect to the objective function f : X → R to analyze in the global space X ,
2. γ-set X γ : A set of configurations X γ that achieves the top-γ quantile in the global space X, and
3. Marginal γ-set PDF pd(xd|X γ ): The marginal PDF of the γ-set PDF p(x|X γ ):

We provide the formal definitions in Appendix B and the conceptual visualizations in Figure 2. Note that X := X1 ×· · ·× XD is the search space, Xd ⊆ R for d ∈ [D] := {1, . . . ,D} is the domain of the d-th HP, xs ∼ Xs ⊆ R|s| denotes xs is sampled from the uniform distribution on Xs where s ⊆ [D], and X-d ⊆ RD-1 is X without the d-th dimension. Furthermore, we consistently denote the PDF of the uniform distribution as uniform PDF and follow the assumptions stated in Appendix C.1.

2.2 f-ANOVA

In this section, we describe f-ANOVA for one-dimensional effects and refer to more details about the general version in Appendix A.3. Suppose we would like to quantify HPI of a function f(x) defined on X, then global HPI [Hooker, 2007] requires (see Figure 3 for the intuition):

1. Global mean:

2. Marginal mean:

(b) The conceptual visualization of global HPI on f (x1, x2) = x21 + x22/10.

Figure 3: The conceptual visualization of global HPI on f(x1, x2) = x21 + x22/10. Left: the landscape of f(x1, x2) (darker is better) and the marginal means f1(x1) and f2(x2). The color does not change a lot in the vertical direction (the x2-axis) while it does in the horizontal direction (the x1-axis). Right: the marginal means and the mean value in the same plane. As the marginal mean f1 (the red dashed line) has a large variance, x1 is more important. On the other hand, as the marginal mean f2 (the blue dashed line) has a small variance, x2 is less important.

3. Marginal variance:

Note that f(x|xd) implies that we fix the d-th HP of x to xd. When we denote the global variance as v0, the ratio vd/v0 is the global HPI of the d-th HP and in essence, the magnitude of the marginal variance represents the relative importance.

2.3 Local f-ANOVA in PriorWorks

To the best of our knowledge, there are two papers that mention local HPI (and both use f-ANOVA). Hutter et al. [2014] mentioned local HPI can be quantified by taking:

however, since this measure is biased depending on the global space design as discussed in Appendix B.4.1, we need to consider the integral only over a local space as stated in Section 3.1. Biedenkapp et al. [2018] proposed the following

where xopt ∈ RD is the optimized setting and xopt -d ∈ RD-1 is xopt without the d-th dimension. The authors mention that this measure is a local HPI measure; however, this measure is also not a local HPI measure in our definition and we show that this is the case using a toy example in Appendix B.4.2.

3 Local f-ANOVAUsing Pearson Divergence

In this section, we first provide the definition of local HPI and describe how to define a local space. For simplicity, we name this local space definition as Lebesgue split 1 . Then we introduce fast algorithm using PED between two KDEs to compute local HPI and benchmark the speed of the algorithm compared to the f-ANOVA implementation based on random forests [Hutter et al., 2014]. Notice that since higher orders of HPI require exponential amounts of computations and usually lack interpretability, our discussion does not focus on higher orders; however, we derive the formula for higher orders and show them in Eqs. (47), (52) in Appendix C. The theoretical details for this section are available in Appendix B.

1 The name comes from the fact that we define a local space by a function value as in the Lebesgue integral in contrast to the definition of a local space by bounds for each dimension, which we name Riemann split.

3.1 Local Hyperparameter Importance

In this section, we assume that we have a set of (sorted) observations D := {(xn, f(xn))}Nn=1 such that f(x1) ≤ · · · ≤ f(xN). Then, the top-γ-quantile observations are Dγ = {(xn, f(xn))}⌈γN⌉ n=1 and the γ-quantile value is f γ := f(x⌈γN⌉).

Local Space Defined by Lebesgue Split

To begin with, we formally define local HPI:

Definition 1 (Local HPI) Given a subspace X ⋆ ⊆ X, local HPI is global HPI vd/v0 in the subspace X ⋆ .

Recall that vd is the marginal variance of the d-th HP and v0 is the global variance. Based on Definition 1, the prior works on making f-ANOVA local discussed in Section 2.3 are not local HPI measures; see Appendix B.4 for more details. As our local HPI obviously depends on the choice of X ⋆ , local HPI is a very general concept; therefore, we focus on the socalled Lebesgue split to specify a local space in this paper. In the Lebesgue split, we obtain a local space X ⋆ as follows:

1. Fix a threshold f ⋆ (we use the γ-quantile value f γ in this paper instead), and
2. Obtain the sublevel set X ⋆ := {x ∈ X|f(x) ≤ f ⋆ } based on f ⋆ (X ⋆ becomes the γ-set X γ when f ⋆ = f γ ).

Recall that the definitions of the γ-quantile value and the γset are available in Section 2.1. More intuitively, the red domains in Figure 2 are the local space of each example. In this paper, we use f(x⌈γN⌉) as f γ . The advantages of the Lebesgue split are to:

1. require only one parameter f ⋆ while the Riemann split, which we split along each dimension by specifying bounds, requires at least 2×D parameters,
2. be able to focus on the analysis in promising domains where we are interested, and
3. be able to remove the sampling bias caused by a nonuniform sampler when using the formula of local HPI.

We further discuss the strengths and drawbacks of the Lebesgue split compared to the Riemann split in Appendix B.5.

Formula of Local Hyperparameter Importance

Now we discuss the computation of local HPI. In Eqs. (2)- (4), we take the expectation over the uniform distribution of the global space X. In the same vein, it is natural to consider the expectation over the uniform distribution of the local space X γ for local HPI as well. Although the computation is not obvious, we can compute the expectation of a measurable function f(x) over the local space X γ if we use the following trick:

(b) The conceptual visualization of global HPI by PED.

Figure 4: The conceptual visualization of global HPI using PED on f(x1, x2) = x21 + x22/100. As we consider global HPI in this example, the γ-set PDF is the uniform PDF. Left: the landscape of f(x1, x2) (darker is better) and the marginal γ ′ -set PDFs p1(x1|Dγ′ ) (= p1, the red dashed line) and p2(x2|Dγ′ ) (= p2, the blue dashed line). The dots represent observations (datasets) and the red dots achieve the top-γ ′ quantile. p1, p2 are estimated by Eq. (15) using the red dots. Right: the marginal γ ′ -set PDFs in the same plane. While p1 (the red dashed line) sharply peaks at the center, p2 (blue dashed line) is close to the uniform PDF. This implies that the latter is closer to the uniform PDF, and thus x2 is less important.

where γ = Ex∼X [I[x ∈ X γ ]] is a normalization constant. Recall that x ∈ X γ and f(x) ≤ f γ are equivalent. Similarly, the marginal mean of f(x) is computed as follows:

where V γ d (xd) := Ex-d∼X-d [I[x ∈ X γ |xd]] is a normalization constant. Then local HPI is generally computed as follows:

1. Local mean:

2. Local marginal mean:

3. Local marginal variance:

Note that xd ∼ V γ d implies that xd is sampled from the distribution of the PDF V γ d (xd)/Z where Z ∈ R+ is the normalization constant. As the series of computations requires a Monte-Carlo sampling in a D - 1 dimensional space, the time complexity incurs the curse of dimensionality. In the next section, we introduce fast algorithm to compute local HPI in exchange for the scale ignorance.

Algorithm 1 Local PED-ANOVA

• D = {(xn, f(xn))}Nn=1 (Dataset to analyze), γ, γ ′ (User-defined quantiles of top domains)
• 1: ▷ See Appendices E.2, E.3 for practical usages
• 2: Sort D in ascending order by f
• 3: ▷ |Dγ | ≥ 2 and |Dγ ′ | ≥ 2 must hold
• 4: Pick the top-γ and -γ ′ quantile observations Dγ ,Dγ ′
• 5: for d = 1, . . . ,D do
• 6: Count occurrences of unique values c(n) d
• 7: Build KDEs pd(·|Dγ ), pd(·|Dγ ′ ) by Eq. (15)
• 8: Compute vγ d by Eq. (16)
• 9: return {vγ d }Dd=1

3.2 Fast Algorithm by Pearson Divergence

If we analyze the binary function b(x|X γ ′ ) := I[x ∈ X γ ′ ] = I[f(x) ≤ f γ ′ ] instead of f(x), HPI can be efficiently computed where γ ′ (< γ) is another quantile to define the binary function in the local space X γ . First, we prove the following theorem:

Theorem 1 Given the binary function b(x|X γ ′ ) and the γ ′ - and γ-set PDFs p(x|X γ ′ ), p(x|X γ ) where γ ′ < γ, the local marginal variance of each dimension d ∈ [D] is:

The proof is provided in Appendix C.3 and higher orders of HPI can be computed by Eq. (52) in Appendix C.3. Note that PED between the PDFs p, q defined on Xd is computed as:

As we do not have the ground truth of the marginal γ ′ - and γ-set PDFs, we replace them with KDEs. The tricks of this computation are that (1) the marginal γ-set PDF can be easily estimated by (1D) KDE as follows and (2) we only need to take the average in 1D space:

Note that xn,d ∈ Xd is the d-th dimension of xn and k is a kernel function. Although the query of this function still requires O(N), the time complexity scales down to O(nd) where nd ∈ Z+ is the number of unique values in the d-th HP if we use the following compression:

where x(n) d is the n-th unique value in the d-th HP and c(n) d is the occurrences of this value in Dγ . Note that we discretize a continuous HP xd ∈ [L,R](L < R) (if exists) as xd ∈ {L + n(R - L)/(nd - 1)}nd-1 n=0 to apply Eq. (15)

and the discretization error of marginal variances is bounded by O( 1 nd ) under some assumptions. Since we can avoid Monte-Carlo samplings with the discretization and the total time complexity is reduced to O(N +n2d ), this is a trade-off; see Proposition 2 in Appendix C.4 for more details. Hence Eq. (11) is approximated as the following closed-form:

where Z := ∑nd n=1 pd(x(n) d |Dγ ) is a normalization constant and the time complexity of this computation is O(n2d ). Algorithm 1 shows the pseudocode for the local HPI computation. Note that global HPI, whose computation is detailed in Appendix B.2, can be computed by replacing pd(·|Dγ ) with the uniform PDF in the d-th dimension ud(xd) (or pd(·|Dγ ) with γ = 1, i.e. pd(·|D), as discussed in Appendix E.2 when collecting D by a non-uniform sampler). Figure 4 presents an example of global HPI with our method on a 2D toy function.

4 Performance Validation

4.1 Setup

In this section, we consistently use the following function:

where xd ∈ [-5, 5] for all d ∈ {1, 2, 3, 4} and the weights wd : R→W follow:

and W := {Wd}3d=0 = {50 , 5-1 , 5-2 , 5-3 }. This function has different trends of HPI in global and local spaces. While the order of HPI is x1, x2, x3, x4 in the global space, it is x2, x3, x4, x1 in the local space {x ∈ X |∀d ∈ [4], |xd| < 1}.

In this experiment, we discretized the HPs with nd = 1001 and all samples were drawn from the uniform distribution. Furthermore, all experiments were run on the hardware with Intel Core i7-10700 and we used the f-ANOVA implementation with the default parameter setting by Optuna 2 . Note that the Optuna implementation is based on Hutter et al. [2014].

4.2 Effect of Scale Ignorance in Global HPI

Since PED-ANOVA uses I[f(x) ≤ f ⋆ ] instead of f(x), it cannot capture scale information. On the other hand, since our objective is to identify important HPs, we would like to test if PED-ANOVA can identify important HPs. In the experiment, we used γ ′ = 0.1. Figure 5 shows the cumulative global HPI ratio of each method. As seen in the figure, while both methods could identify the most important HP x1, we can see the difference in the HPI of x2. While PED-ANOVA tells us x2 has about 20% of contribution to achieve the top10%, f-ANOVA tells us x2 has about 3% of contribution. This difference comes from whether we ignore the scale of the objective function or not. Since f-ANOVA considers scale and it magnifies the contribution in the tail of the function, it dilutes the HPI of x2, which has less weight in the tail. Note that "tail" refers to the domains that cause a lot of variations, yet not critical for the final result, in the objective function f and |xd| ≥ 1 is the tail in our case; more details in Appendix B.6. On the other hand, the HPI of x1 by PED-ANOVA is not strongly biased by the tail due to the scale ignorance nature. This leads to more importance in x2. Although our method loses scale information, the ignorance of scale allows us to abandon the information from the tail and focus only on the information from the promising domain, which is γ ′ - set in our case. Furthermore, this remarkable property makes the meaning of HPI, which is how important each HP is to achieve the top-γ ′ quantile, very clear and practitioners can extract the nuance of each HP for specific local spaces.

4.3 Query Speed

As mentioned previously, one of the benefits of our method is the query speed, and we would like to benchmark how quick our method is in this section. In the experiment, we used γ ′ = 0.1. Figure 6 presents the runtime with respect to the number of data points. While f-ANOVA requires more than a minute with 103 data points and more than a week with 105 data points, our method provides the results in a minute even with 108 data points. In Appendix D.2, we test our method with higher dimensionality to see the number of data points required for convergence.

2 https://github.com/optuna/optuna

4.4 Local Importance Measure

Finally, we check if our method can successfully identify important HPs in promising domains. The objective function f(x) is designed so that while x1 is important and x3 is trivial in the global space, x1 is less important and x3 is important in the local space. The goal of this experiment is to check whether our method can provide this insight. In the experiment, we used N = 104 .

Figure 7 shows the results. As discussed already, global HPI could identify the order of HPI appropriately. For local HPI, our method could tell us that x2, x3 are the most important HPs in the local space and x1 is less important as expected. Note that since the γ = 0.1-set already narrows down the domain of x2, but not x3, this dilutes the HPI of x2 and increases the HPI of x3. Prior works cannot provide this interpretation as discussed in Appendix B.4.

5 Real-World Usecase by JAHS-Bench-201

5.1 Setup

In order to further verify our proposed algorithm against a real-world application, we applied PED-ANOVA to analyze the search space of JAHS-Bench-201 [Bansal et al., 2022], which is a surrogate benchmark for HPO and has a very large search space in the context of extant HPO benchmarks. We constructed the dataset D in Algorithm 1 by querying JAHSBench-201 for the validation accuracy, i.e. f(x), of N lattice points, where N = 41,343,750, generated by discretizing the JAHS-Bench-201 search space (see Table 2 of Appendix D.1). Although JAHS-Bench-201 can be queried for model performance metrics on 3 different image classification datasets, for the sake of brevity, we discuss only the experiments performed on CIFAR10 here and include the results on the other datasets in Appendix D.3. Due to the computational complexity of f-ANOVA, we could use only 104

HPI ratio (%)

data points for it, in contrast to PED-ANOVA, and calculated the mean of HPI over 10 independent runs. Since the surrogate models in JAHS-Bench-201 are trained XGBoost models and XGBoost's outputs are deterministic, we query each lattice point only once. In the analysis, we would like to answer the following research questions (RQs):

• RQ1: Does global HPI of our method provide the same important HPs as f-ANOVA with Eq. (5)?
• RQ2: Is the scale ignorance necessary for matching the intuition of achieving the top-γ ′ quantile?
• RQ3: Does local HPI help detect potentially important HPs or trivial HPs?

In order to answer RQs, we provide Table 1 with the HPI of each HP in CIFAR10 of JAHS-Bench-201 and Figure 8 to visualize the γ = 0.1- and γ = 0.01-set PDFs as the blue and the red shadows, respectively. Strictly speaking, discrete probability distributions are not PDFs due to discrete space; however, we use the term γ-set PDF for the sake of consistency. We applied f-ANOVA to f(x) (Normal), min(f(x), f γ ′=0.1 ) (Global 0.1), and min(f(x), f γ ′=0.01 ) (Global 0.01), and PED-ANOVA

with γ ′ = 0.1 (Global 0.1), γ ′ = 0.01 (Global 0.01), and γ = 0.1, γ ′ = 0.01 (Local). Recall that Global 0.1 and Global 0.01 for f-ANOVA are based on Eq. (5) and f γ ′=0.1 means f⌈|D|/10⌉ given a dataset sorted by fn. Although we used the uniform PDF to compute global HPI in this experiment, practitoners should use pd(·|D)(γ = 1) instead of the uniform PDF for the post-hoc analysis of HPO when using a non-uniform sampler (e.g. Bayesian optimization) to remove sampling bias as discussed in Appendix E.2.

5.2 Analysis & Interpretation

To answer RQ1, we compare the column (Global 0.1, Ours) to (Global 0.1, Original) and the column (Global 0.01, Ours) to (Global 0.01, Original) in Table 1. We observe that both PED-ANOVA and f-ANOVA indicated the same top-2 important HPs although the 3rd -best HPs slightly varied. This result further verifies the validation in Section 4.2.

To answer RQ2, we discuss the results of (Global 0.1, Ours) and (Global 0.01, Ours) in the context of (Normal, Original) to assess the impact that the tail of f(x) discussed in Section 4.2 has on f-ANOVA. The most important takeaway from this comparison is the misclassification of Op.3 as the most important and of TrivialAugment as the least important HP to optimize over by the original f-ANOVA in the global setting. As can be verified by looking at the γ-set PDFs, even for γ = 0.01, Op.3's values are distributed very evenly even when TrivialAugment and Width multiplier have already shown convergence. This clearly indicates that Op.3 is not very important to optimize for achieving the top-1% performance and may or may not become relevant in even higher quantile regimes. At the same time, both columns' values agree on the importance of Op.1 and Op.6. Therefore, to answer RQ2, scale invariance indeed helps to successfully identify HPI for HPs that would have been misclassified by the (Normal, Original) setting.

Finally, for RQ3, we compare the column (Global 0.01, Ours) to (Local, Ours). We observe that the HPI of Width multiplier drops sharply from the Global 0.01 setting to the Local setting. Simultaneously, the HPI of TrivialAugment increases sharply across the same. This suggests that optimizing Width multiplier is no longer important when moving from the top-10% to the top-1% performance but optimizing TrivialAugment is very important. The reason behind this change becomes clear when we observe the change in γ-set PDFs of the two HPs in Figure 8. Both the γset PDFs for Width multiplier are sharply peaked at 16, indicating that no further optimization is needed on Width multiplier. However, the γ-set PDFs for TrivialAugment only start peaking at the value True for the γ = 0.01-set PDF. This clearly demonstrates that local HPI is necessary for deriving the correct interpretation in the top-γ ′ quantiles, since (Local, Ours) successfully identifies the relative importance of optimizing the two HPs. Last but not least, if both global and local HPI with wished quantiles γ, γ ′ exhibits low values, removing such HPs, e.g. Activation function, is expected to have a less negative impact although it is insecure to remove HPs, e.g. Op.1, only by looking at global HPI.

6 Conclusions

In this paper, we reformulated f-ANOVA for local HPI and introduced the fast algorithm to compute local HPI by PED. In the series of experiments on a toy function, we confirmed that our method can quantify both global and local HPI appropriately, and efficiently compute HPI in a second with 105 data points while the prior work takes several days. In the analysis of JAHS-Bench-201, we provided a concrete example of how to use our method on benchmark datasets and showed that only using global HPI could be misleading. Due to the space limit, we defer a discussion of practical usecases and limitations of our method to Appendix E. Our implementation is available at https://github.com/nabenabe0928/local-anova/.

Acknowledgments

The authors appreciate the valuable contributions of the anonymous reviewers. Robert Bosch GmbH is acknowledged for financial support. The authors also acknowledge funding by European Research Council (ERC) Consolidator Grant "Deep Learning 2.0" (grant no. 101045765). Views and opinions expressed are however those of the authors only and do not necessarily reflect those of the European Union or the ERC. Neither the European Union nor the ERC can be held responsible for them.