{'fraction_non_alphanumeric': 0.03272971723877148, 'fraction_numerical': 0.008802151170118903, 'mean_word_length': 4.6455170777988615, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Active Assessment of Prediction Services as Accuracy Surface Over Attribute Combinations

Our goal is to evaluate the accuracy of a black-box classification model, not as a single aggregate on a given test data distribution, but as a surface over a large number of combinations of attributes characterizing multiple test data distributions. Such attributed accuracy measures become important as machine learning models get deployed as a service, where the training data distribution is hidden from clients, and different clients may be interested in diverse regions of the data distribution. We present Attributed Accuracy Assay (AAA)--a Gaussian Process (GP)--based probabilistic estimator for such an accuracy surface. Each attribute combination, called an 'arm', is associated with a Beta density from which the service's accuracy is sampled. We expect the GP to smooth the parameters of the Beta density over related arms to mitigate sparsity. We show that obvious application of GPs cannot address the challenge of heteroscedastic uncertainty over a huge attribute space that is sparsely and unevenly populated. In response, we present two enhancements: pooling sparse observations, and regularizing the scale parameter of the Beta densities. After introducing these innovations, we establish the effectiveness of AAA in terms of both its estimation accuracy and exploration efficiency, through extensive experiments and analysis.


Introduction
Increasing concentration of big data and computing resources has resulted in widespread adoption of machine learning as a service (MLaaS). The best-performing NLP, speech, image and video recognition tools are now provided as network services. MLaaS comes with few accuracy specifications or service level agreements, perhaps only leaderboard numbers from benchmarks that may not be closely related to most clients' deployment data distributions. The client, therefore, finds it difficult to choose the best provider without extensive pilot trials [1]. Different clients may need to deploy the service on very different data distributions, with possibly widely different accuracy.
In such circumstances, we propose that a service provider, or a service standardization agency, publish the accuracy of the classifier, not as one or few aggregate numbers, but as a surface defined on a space of input instance attributes that capture the variability of consumer expectations. Indoor/outdoor, day/night, urban/rural may be attributes of input images for visual object recognition tasks. Speaker age, gender, ethnicity/accent may be attributes of input audio for speech recognition tasks. We call a combination of attributes in their Cartesian space an arm (borrowing from bandit terminology) 2 . The labeled instances used by the service provider may not represent or cover well the space of attributes of interest to subscribers. Labeled data may be proprietary and inaccessible to prospective consumers and standardization agencies. Whoever estimates the accuracy surface, therefore, needs to actively select instances from an unlabeled pool for labeling, presumably within a restricted budget, to adequately cover the attribute space.
Several recent studies have highlighted the variability in accuracy across data sub-populations [2,3], specifically in the context of fairness [4,5,6], and also proposed active estimation techniques of sub-population accuracy [7,8]. We solve a more general problem where the space of arms (subpopulation) defined by the Cartesian space of attributes grows combinatorially. This inevitably leads to extreme sparsity of labeled instances for many arms. A central challenge is how to smooth the estimate across related arms while faithfully representing the uncertainty for active exploration.
We present Attributed Accuracy Assay (AAA) -a practical system that estimates accuracy, together with the uncertainty of the estimate, as a function of the attribute space. AAA uses these estimates to drive the sampling policy for each attribute combination. Gaussian Process (GP) regression is a natural choice to obtain smooth probabilistic accuracy estimates over arm attributes. However, a straightforward GP model fails to address the challenge of heteroscedasticity that we face with uneven and sparse supervision across arms. We model arm-specific service accuracy as drawn from a Beta density that is characterized by mean and scale parameters, which are sampled from two GPs that are informed by suitable trained kernels over the attribute space. We propose two further enhancements to the training of this model. First, we recognize an over-smoothing problem with GP's estimation of the Beta scale parameters, and propose a Dirichlet likelihood to supervise the relative values of scale across arms. Second, we recognize that arms with very low support interfere with learning the kernel parameters of the GPs. We mitigate this by pooling observations across related arms. With these fixes, AAA achieves the best estimation performance among competitive alternatives.
Another practical challenge in our setting is that some attributes of instances are not known exactly. For example, attributes, such as camera shutter speed or speaker gender, may be explicitly provided as meta information attached with instances. But other attributes, such as indoor/outdoor, or speaker age, may have to be estimated noisily via another (attribute) classifier, because accurate human-based acquisition of attributes would be burdensome. AAA also tackles uncertain attribute inference. Its attribute classifiers are trained on a small amount of labeled data and their error rates are modeled in a probabilistic framework.
We report on extensive experiments using several real data sets. Comparison with several estimators based on Bernoulli arm parameters, Beta densities per arm, and even simpler forms of GPs on the arm Beta distributions, shows that AAA is superior at quickly cutting down arm accuracy uncertainty.
Summarizing, our contributions are: • We motivate and define the problem of accuracy surface estimation over a large space of attribute combinations. • Our proposed estimator AAA fits a Beta density for every attribute combination (arm), with its parameters smoothed via two GPs to capture heteroscedastic uncertainty of each arm's accuracy under limited data settings. • We propose two important components included in AAA: 1) a Dirichlet regularization to control over-smoothing of the Beta scale parameters, and 2) pooled observations to reduce over-fitting of a GP-associated kernel to sparse arms. • We show significant gains in terms of both estimation quality and the efficiency of exploration on four real classification models compared to existing methods. AAA obtains an average 80% reduction in macro averaged square error over the existing methods.

Problem Setup
Our goal is to evaluate a given machine learning service model S used by a diverse set of consumers. The service S : X → Y could be any predictive model that, for an input instance x ∈ X , assigns an output labelŷ ∈ Y, where Y is a discrete label space. Let y(x) denote the true label of x and Agree(y,ŷ) denote the match between the two labels. For scalar classification, Agree(y,ŷ) is in {0,1}. For structured outputs, e.g., sequences, we could use measures like BLEU scores in [0,1]. Classifiers are routinely evaluated on their expected accuracy on a data distribution P (X , Y): (1) We propose to go beyond this single measure and define accuracy as a surface over a space of attributes of the input instances. Let A denote a list of K attributes that capture the variability of consumer expectation on which the service S will be deployed. For instance, visual object recognition is affected by the background scene, and facial recognition is affected by demographic attributes. We use A(x) ∈ A to denote the vector of values of attributes of input x and A to denote the Cartesian product of the domains of all attributes. An attribute could be discrete, e.g., the ethnicity of a speaker; Boolean, e.g., whether a scene is outdoors/indoors; or continuous, e.g., the age of the speaker in speech recognition. Some of the attributes of x, for example the camera settings of an image, may be known exactly, and others may only be available as a distribution M k (a k |x) for an attribute a k ∈ A, obtained from a pre-trained probabilistic classifier.
Generalizing from a single global expected accuracy (1), we define the accuracy surface ρ : A → [0, 1] of a service S at each attribute combination a ∈ A, given a data distribution P (X , Y), as Our goal is to provide an estimate of ρ(a) given two kind of data sampled from P (X , Y): a small labeled sample D, and a large unlabeled sample U . In addition, we are given a budget of B instances for which we can seek labels y from a human by selecting them from U . Applying M k to all of U is, however, free of cost.
We aim to design a probabilistic estimator for ρ(a), which we denote as P (ρ|a) where ρ ∈ [0, 1] and a ∈ A. This is distinct from active learning, which selects instances to train the learner toward greater accuracy, and also active accuracy estimation [7], which does not involve a surface over as.
We also show that standard tools to regress from a to ρ are worse than our proposal.
We measure the quality of our estimate as the square error between the gold accuracy ρ(a) and the mean of the estimated accuracy distribution P (ρ|a). Our estimator distribution naturally gives an idea of the posterior variance of accuracy estimate of each attribute combination, which we use for uncertainty-based exploration.

Proposed Estimator
We will first review recent work that leads to candidate solutions to our problem, discuss their limitations, and finally present our solution. Initially, to keep the treatment simple, we assume A(x) and gold y (hence c = Agree(S(x), y), the service correctness bit) is known for all instances. Later in this section, we remove these assumptions.
The simplest option is to ignore any relationship between arms, and, for each arm a, fit a suitable density over ρ(a). When this density is sampled, we get a number in [0, 1], which is like a coin head probability used to sample correctness bits c. For representing uncertainty of accuracy values (which are ratios between two counts), the Beta distribution B(·, ·) is a natural choice. We call this baseline method Beta-I.
The variance of the estimated Beta density can be used for actively sampling arms. Ji et al. [7] describe a related scenario, stressing on active sampling. However, this approach cannot share observations or smooth the estimated density at a sparsely-populated arm with information from similar arms. In our real-life scenario, we expect accuracy surface smoother and the number of arms to be large enough that many arms will get very few, if any, instances.
The second baseline method, which we call BernGP, is to view the (a, c) instances in D as a standard classification data set with the binary c values as class label and a as input features. Given the limited data, we can use the well-known GP classification approach [9] for fitting smooth values ρ as a function of a. Suppose the arms a can be embedded to V(a) in a suitable space induced by some similarity kernel. In this embedding space, we expect the accuracy of S to vary smoothly. Given a kernel K 1 (a, a ) to guide the extent of sharing of information across arms, a standard form of this GP would be The GP can give estimates of uncertainty of ρ(a), which may be used for active sampling of arms.
As we will demonstrate, such GP-imposed estimate of uncertainty of ρ(a) is inadequate, because it loses sight of the number of supporting observations at each arm, which could be very diverse. This is because the standard GP assumption of homoscedasticity, that is, identical noise around each arm is violated when observations per arm differ significantly. We therefore need a mechanism to separately account for the uncertainty at each arm, even the unexplored ones, to guide the strategy for actively collecting more labeled data.

The basic BetaGP proposal
We model arm-specific noise by allowing each arm to represent the uncertainty of ρ a , not just by an underlying GP as in BernGP above, but also by a separate scale parameter. Further, the scale parameter is smoothed over neighboring arms using another GP. The influence of this scale on the uncertainty of ρ a is expressed by a Beta distribution as follows: where we use φ(•), ψ(•) to denote the parameters of the Beta distribution at arm a. The Beta distribution is commonly represented via α, β parameters whereas we chose the less popular mean (φ) and scale (ψ) parameters. While these two forms are functionally equivalent with φ = α α+β , ψ = α + β, we preferred the second form because imposing GP smoothness across arms on the mean accuracy and scale seemed more meaningful. We validate this empirically in the Appendix B.
Two kernel functions K 1 (a, a ), K 2 (a, a ) defined over pairs of arms a, a ∈ A control the degree of smoothness among the Beta parameters across the arms. We use an RBF kernel defined over learned shared embeddings V(a): where s 1 , s 2 , l 1 , l 2 denote the scale and length parameters of the two kernels. The scale and length parameters are learned along with the parameters of embeddings V(a) during training.
Initially, we assume we are given a labeled dataset D = {(x i , a i , y i ) : i = 1 . . . , I} with attribute information available. Using predictions from the classification service S, we associate a 0/1 accuracy c i = Agree(y i , S(x i )). We can thus extend D to Let c a = i:A(xi)=a c i denote the total agree score in arm a. Let n a denote the total number of labeled examples in arm a. The likelihood of all observations given functions f, g decomposes as a product of Beta-binomial 3 distributions at each arm as follows: where B is the Beta function, and the second expression is a rewrite of the Beta-binomial likelihood.
During training we calculate the posterior distribution of functions f, g using the above data likelihood Pr(D|f, g) and GP priors given in eqns. (5) and (6). The posterior cannot be computed analytically given our likelihood, so we use variational methods. Further, we reduce the O(|A| 3 ) complexity of posterior computation, using the inducing point method of Hensman et al. [9], whereby we learn m locations u ∈ R d×m , mean µ ∈ R m , and covariance Σ ∈ R m×m of inducing points. Doing so brings down the complexity to O(m 2 |A|), m |A|. These parameters are learned end to end with the parameters of the neural network used to extract embeddings V(a) of arms a, and kernel parameters s 1 , s 2 , 1 , 2 . We used off-the-shelf Gaussian process library: GPyTorch [10] to train the above likelihood with variational methods. Details of this procedure can be found in the Appendix C. We denote the posterior functions as P (f |D), P (g|D). Thereafter, the mean estimated accuracy for an arm a is computed as We call this setup BetaGP. Next, we will argue why BetaGP still has serious limitations, and offer mitigation measures.

Supervision for scale parameters
We had introduced the second GP g a to model arm-specific noise, and similar techniques have been proposed earlier by Lázaro-Gredilla and Titsias [11], Kersting et al. [12], Goldberg et al. [13], but for heteroscedasticity in Gaussian observations. However, we found the posterior distribution of scale values ψ(g a ) at each arm tended to converge to similar values, even across arms with orders of magnitude difference in number of observations n a . On hindsight, that was to be expected, because the data likelihood (8) increases monotonically with scale ψ a . The only control over its converging to ∞ is the GP prior g ∼ GP (0, K 2 ). In Appendix D, we illustrate this phenomenon with an example. We propose a simple fix to the scale supervision problem. We expect the relative values of scale across arms to reflect the distribution of the proportion of observations na n across arms (with n = a n a ). We impose a joint Dirichlet distribution using the scale of arms ψ(g a ) as parameters, and write the likelihood of the observed proportions as (with Γ denoting Gamma function): We call this BetaGP-SL. With this as an additional term in the data likelihood, we obtained significantly improved uncertainty estimates at each arm, as we will show in the experiment section.

Pooling for sparse observations
Recall that the observations are accumulation of 1/0 agreement scores for all instances that belong to an arm. Given the nature of our problem, arms have varying levels of supervision, and also highly varying true accuracy values. Even when the available labeled data is large, many arms will continue to have sparse supervision because they represent rare attribute combinations. The combination of high variance observations and sparse supervision could lead to learning of non-smooth kernel parameters. In Appendix D, we demonstrate with a simple setting that GP parameters learned on noisy observations under-represent the smoothness of the surface. The situation is further aggravated when learning a deep kernel. This problem has resemblance to "collapsing variance problem" [14] such as when Gaussian mixture models overfit on outliers or when topic models overfit a noisy document in the corpus. Instead of depending purely on GP priors to smooth over these noisy observations, we found it helpful to also externally smooth noisy observations. For each arm a with observations below a threshold, we mean-pool observations from some number of nearest neighbors, weighted by their kernel similarity with a. We will see that such external smoothing resulted in significantly more accurate estimates particularly for arms with extreme accuracy values. We call this method BetaGP-SLP (note that this also includes the scale supervision objective described in the previous section). Two other mechanisms take us to the full form of the AAA system, which we describe next.

Exploration
The variance estimate of an arm informs its uncertainty and is commonly used for efficient exploration [15]. Let P (f |D), P (g|D) denote the learned posterior distribution of the GPs. Using these, the estimated variance at an arm is given as: where the expected value is given in eqn. (10). We use sampling to estimate the above expectation. The arm to be sampled next is chosen as the one with the highest variance among unexplored arms. We then sample an unexplored example with highest affiliation (P (a | x)) with the chosen arm.

Modeling Attribute Uncertainty
Recall that attributes of an instance x are obtained from models M k (a k |x), k ∈ [K], which may be highly noisy for some attributes. Thus, we cannot assume a fixed attribute vector A(x) for an instance x. We address this by designing a model that can combine these noisy estimates into a joint distribution P (a|x) using which, we can fractionally assign each instance x i across arms. A baseline model for P (a|x) would be just the product K k=1 M k (a k |x). However, we expect values of attributes to be correlated (e.g. attribute 'high-pitch' is likely to be correlated with gender 'female'). Also, the probabilities M k (a k |x) may not be well-calibrated.
We therefore propose an alternative joint model that can both recalibrate individual classifiers via temperature scaling [16], and model their correlation. We have a small seed labeled dataset D with gold attribute labels, independent noisy distributions from each attribute model M k (a k |x), and an unlabeled dataset U . We prefer simple factorized models. We factorize log Pr(a|x) as a sum of temperature-weighted logits and a joint (log) potential as shown in expression (13) below.
log Pr(a|x) = log Pr(a 1 , a 2 , · · · , a K |x) = Here N denotes a dense network to model the correlation between attributes, and t 1 , . . . , t K denote the temperature parameters used to rescale noisy attribute probabilities. The maximum likelihood Z i denotes the partition function for an example x i which requires summation over A. Exact computation of Z i could be intractable especially when A is large. In such cases, Z i can be approximated by sampling. In our case, we could get exact estimates.
In addition to D, we use the unlabeled instances U with predictions from attribute predictors filling the role of gold-attributes. Details on how we train the parameters on large but noisy U and small but correct D can be found in the Appendix E.
The estimation method of BetaGP-SLP with variance based exploration and calibration described here constitute our proposed estimator: AAA. Detailed pseudo-code of AAA is given in the Appendix I.

Experiments
Our exploration of various methods and data sets is guided by the following research questions.

Data sets and tasks
We experiment with two real data sets and tasks. Our two tasks are male-female gender classification with two classes and animal classification with 10 classes.
Male-Female classification (MF): CelebA [17] is a popular celebrity faces and attribute data set that identifies the gender of celebrities among 39 other binary attributes. The label is gender. The accuracy surface spans various demographic, style, and personality related attributes. We hand-pick a subset of 12 attributes that include attributes that we deem important for gender classification among some other gender-neural attributes such as if the subject is young or wearing glasses (see Appendix F for more details). We used a random subset of 50,000 examples from the dataset for training classifiers on each of the 12 attributes using a pretrained ResNet-50 model. The remaining 150,000 examples in the data set are set as the unlabeled pool from which we actively explore new examples for human feedback. The twelve binary attributes make up for 2 12 = 4, 096 attribute combinations.
Animal classification (AC): COCO-Stuff [18] provides an image collection. For each image, labels for foreground (cow, camel) and background (sky, snow, water) 'stuff' are available. Visual recognition models often correlate the background scene with the animal label such as camel with deserts and cow with meadows. Thus, foreground labels are our regular y-labels while background stuff labels supply our notion of attributes.
We collapse fine background labels into five coarse labels using the dataset provided label hierarchy. These are: water, ground, sky, structure, furniture. The Coco dataset has around 90 object (foreground)   labels. Here we use a subset of 10 labels corresponding to animals. We take special care to filter out images with multiple/no animals and adapt the pixel segmentation/classification task to object classification (see the Appendix F for more details). The image is further annotated with the five binary labels corresponding to five coarse stuff labels. The scene descriptive five binary labels and ten object labels make up for 32×10 = 320 attribute combinations.

Service Models
For the MF task, we use two service models (S). MF-CelebA is a service model for gender classification. To simulate separate D and U , it is trained on a random subset of CelebA with a ResNet50 model. MF-IMDB is a publicly available 4 classifier trained on IMBD-Wiki dataset, also using the ResNet50 architecture. The attribute predictors are trained using ResNet50 on a subset of the CelebA dataset for both service models.
For the AC task, we use two publicly available 5 service models (S). AC-COCOS was trained on COCOS data set with 164K examples. AC-COCOS10k was trained on COCOS10K, an earlier version of COCOS with only 10K instances. We use these architectures for both label and attribute prediction. See Appendix F, G for more details on accuracy surface, attribute predictor, service models and their architecture. In Figure 1, we illustrate some statistics of the shape of the accuracy surface for the four dataset-task combinations. Although S's mean accuracy (rightmost bars) is reasonably high, the accuracy of the arms in the 10% quantile is abysmally low, while arms in the top quantiles have near perfect accuracy. This further motivates the need for an accuracy surface instead of single accuracy estimate.

Methods Compared
We compare the proposed estimation method AAA against natural baselines, alternatives, and ablations. Some of the methods, such as Beta-I, BernGP and BetaGP, we have already defined in Section 3. We train methods BernGP and BetaGP using the default arm-level likelihood. We also separately evaluate the impact of our fixes on BetaGP with only scale supervision: BetaGP-SL and along with mean pooling: BetaGP-SLP. We also include a trivial baseline: CPredictor which fits all the arms with a global accuracy estimated using gold D. We do not try sparse observation pooling with Beta-I since there is no notion of per-arm closeness. We also skip it on BernGP since it is worse than BetaGP as we will show below. Recall that Beta-I modeling is related to Ji et al. [7]. Unless otherwise specified, we give equal importance to each arm and report MSE macroaveraged over all arms. Along with macro MSE, we also sometimes report MSE on the subset of 50 worst (true-)accuracy arms, referred to as worst MSE. We report other aggregate errors in the Appendix H. All the numbers reported here are averaged over three runs each with different random seed. The initial set of warm-start examples (D) is also changed between the runs. In the case of BetaGP-SLP, for any arm with observation count below 5, we mean pool from its three closest neighbours.
In the following Sections: 4.5 and 4.6, we compare various estimation and exploration strategies with P (a|x) noise calibrated as described in Section 3.5. In Section 4.7, we study different forms of calibration and demonstrate the superiority of our proposed calibration technique of Equation (13).

Accuracy Estimation Quality
We evaluate methods on their estimation quality when each method is provided with exactly the same (randomly chosen) labeled set. We compare the four service models when fitted on labeled data of size 2,000 and the results appear in Table 2. Note that we only have label supervision on Y in the labeled data. Table 2 shows macro and worst MSE, standard deviation for each metric can be found in Appendix H. In Figure 3, we show worst MSE for a range of labeled data sizes along with their error bars. We make the following observations. Smoothing helps: Since we have a large number of arms, we expect Beta-I to fare worse than its smooth counterparts (BernGP and BetaGP), especially on the worst arms. This is confirmed in the table. In three out of four cases, Beta-I method is worse than even the constant predictor CPredictor on both metrics. Modeling arm specific noise helps: BetaGP is better than BernGP on almost all the cases in the table. Significant gains when the scale supervision problem of BetaGP is fixed: BetaGP-SL is significantly better than BetaGP in the table and figure. Our pooling strategy helps: BetaGP-SLP improves BetaGP-SL over worst MSE without hurting macro MSE as seen in the table and figure.

Exploration Efficiency
We compare different methods that use their own estimated variance for exploring instances to label (Section 3.4), as a function of the number of explored examples -see Figure 4. In most cases, BetaGP-SLP gives the smallest macro MSE, beating Beta-I and BetaGP. Note Beta-I is the exploration method recently suggested in [7]. We observe that BetaGP provides very poor exploration quality, indicating that the uncertainty of arms is not captured well by just using two GPs. In fact, in many cases BetaGP is worse than Beta-I, even though we saw the opposite trend in estimation quality ( Figure 3). These experiments brings out the significant role of Dirichlet scale supervision and pooled observations in enhancing the uncertainty estimates at each arm.

Impact of Calibration
We consider two baselines along with our method explained in Section 3.5: Cal:Raw, which uses the predicted attribute from the attribute models without any calibration and Cal:Temp, which calibrates only the temperature parameters shown in eqn. (13), i.e., without the joint potential part. We refer to our method of calibration using temperature and joint potentials as Cal:Full. We compare these on the four tasks with estimation method set to Beta-I and random exploration strategy. Figure

Related Work
Our problem of actively estimating the accuracy surface of a classifier generalizes the more established problem of estimating a single accuracy score [19,20,21,22,23,24]. For that problem, a known solution is stratified sampling, which partitions data into homogeneous strata and then seeks examples from regions with highest uncertainty and support. If we view each arm as a stratum, our method follows similar strategy. A key difference in our setting is that low support arms cannot be ignored. This makes it imperative to calibrate well the uncertainty under limited and skewed support distribution. The setting of Ji et al. [7] is the closest to ours. However, their work only considers a single attribute which they fit using Beta-I, whereas we focus on the challenges of estimating accuracy over many sparsely populated attribute combinations.
Sub-population performance: Several recent papers have focused on identifying sub-populations with significantly worse accuracy than aggregated accuracy [2,3,6,8,25,26]. Some of these have also proposed sample-efficient techniques [6,8] for estimation of performance on specific sub-groups, such as the ones defined by attributes like gender and race. Our accuracy surface estimation problem can be seen as a generalization where we need to estimate for all sub-groups defined in the Cartesian space of pre-specified semantic attributes. Mitchell et al. [5] recommend reporting model performance under the influence of various relevant demographic/environmental factors as model cards-similar to the accuracy surface.
Experiment design: Another related area is experiment design using active explorations with GPs [27]. Their goal is to find the mode of the surface whereas our goal is to estimate the entire surface. Further, each arm in our setting corresponds to multiple instances, which gives rise to a degree of heteroscedasticity and input-dependent noise that is not modeled in their settings. Lázaro-Gredilla and Titsias [11], Kersting et al. [12] propose to handle heteroscedasticity by using a separate GP to model the variance at each arm. However, we showed the importance of additional terms in our likelihood and observation pooling to reliably represent estimation uncertainty. Wenger et al. [28] propose observation pooling for estimating smooth Betas but they assume a fixed kernel.

Conclusion
We presented AAA, a new approach to estimate the accuracy of a classification service, not as a single number, but as a surface over a space of attributes (arms). AAA models uncertainty with a Beta distribution at each arm and regresses these parameters using two Gaussian Processes to capture smoothness and generalize to unseen arms. We proposed an additional Dirichlet likelihood to mitigate an over-smoothing problem with GP's estimation of Beta distributions' scale parameters. Further, to protect these high-capacity GPs from unreliable accuracy observations at sparsely populated arms, we propose to use an observation pooling strategy. Finally, we show how to handle noisy attribute labels by an efficient joint recalibration method. Evaluation on real-life datasets and classification services show the efficacy of AAA, both in estimation and exploration quality. Limitation and future work: (1) We have evaluated AAA on the order of thousands of arms. Even larger attribute spaces could unearth more challenges.
(2) Identifying relevant attributes for an application can be non-trivial. Future work could devise strategies for attribute selection.
(3) Characterizing test-time data shifts could in itself be hard, particularly for text -there could be subtle changes in word usage, style, or punctuation. A more expressive attribute space needs to be developed for text applications.

Active Assessment of Prediction Services as Accuracy Surface Over Attribute Combinations (Appendix)
Main Section Appendix Source Code Appendix A Section 3.1 Appendix B, C Section 3.2 Appendix D Section 3.5 Appendix E, I Section 4.1, 4.2 Appendix F, G Section 4.4, 4.5 Appendix H Table 6: Mapping between main and appendix sections.

A Source Code
Our code, dataset and instructions for replicating the results can be found at this link.

B Parametric Form of BetaGP
In Section 3.1, we claimed that BetaGP with (mean, scale) parameterization is better than BetaGPαβ with the standard (α, β) parameterization of the Beta distribution. In this section, we present some empirical evidence corroborating the claim.
We compare between the two parametric forms with two service models in Table 7. For BetaGPαβ, we use two GPs, one to approximate the latent value corresponding to α, and other for β. We use soft-plus operation to transform the latent values to their admissible positive α, β values.
We report macro-averaged mean square errors on two tasks in Table 7, when fitting on 2,000 instances -similar to the setting of Section 4.5. We found the BetaGPαβ estimates unstable and far worse, perhaps because smoothness is not expected in either of α, β parameters across arms making the GP's bias ineffective.  Table 7: Comparison of estimation error between BetaGP with (mean, scale) parameterization vs. BetaGPαβ. BetaGPαβ is worse than BetaGP.

C More Details of Gaussian Process (GP) Setup
In this section, we give further details on GP training, posterior approximation and computational cost. This section elaborates on Section 3.1.
In all our proposed estimators, the data likelihood is modeled either by a Bernoulli or a Beta distribution. Data likelihood term of BernGP, BetaGP, are shown in Eqn. (3), Eqn. (9), respectively. Due to the non-Gaussian nature of the data likelihood, the posterior on parameters cannot be expressed in a closed form. Several approximations exist for fitting the posterior especially for the more standard BernGP, we will discuss one such method in what follows. Recall that we model two latent values f, g each modeled by an independent GP. They can be seen to have been drawn from a single GP with even larger dimension and with appropriately defined kernel matrix. For the sake of explanation and with a slight abuse of notation, we denote by f , the concatenation of f and g. The corresponding kernel for the concatenated vector is appropriately made by combining the kernels of either of the latent values with kernel entries corresponding to interaction between f and g set to 0.
Variational methods are popular for dealing with non-Gaussian likelihoods in GP. In this method, we fit a multi-variate Gaussian that closely approximates the posterior, i.e. minimizes D KL (q(f ) P (f |D). GPs are often used in their sparse avatars using inducing points [36,37] that provide approximations to the full covariance matrix with large computation benefits. As a result, q(f) is parameterized by the following trainable parameters (let 'd', 'm' denote the input dimension and number of inducing points resp.): (a) Z, a matrix of size d × m, of locations of 'm' inducing points (b) µ ∈ R m , Σ ∈ R m×m , denoting mean and covariance of the inducing points. In order to minimize D KL (q(f ) P (f |D), a pseudo objective called Evidence Lower Bound (ELBO), shown below, is employed: The first term above in Eqn. (15) maximizes data likelihood, which is Equation 9 in our case. The second term is a regularizer that regresses the posterior fit q(f ) close to the prior distribution P (f ) which is set to standard Normal. We optimize using this objective over all the parameters involved through gradient descent. The required integrals in (15) can be computed using Monte Carlo methods [9]. We describe further implementation details in the next section.

Implementation Details
We use routines from GPytorch 6 library to implement the variational objective. Specifically, we extend ApproximateGP with VariationalStrategy, both of which are GPytorch classes, and set them to learn inducing point locations.
Number of inducing points when set to a very low value could overly smooth the surface and can have high computation overhead when set to a large value. We set the number of inducing points to 50 for all the tasks. The choice of 50 over a larger number is only to ensure reasonable computation speed.
Since we keep getting more observations as we explore, we use the following strategy for scheduling the parameter updates. We start with the examples in the seed set D and update for 1,000 steps. We explore using the variance of the estimated posterior. We pick 12 arms with highest variance and label one example for each arm. After every new batch of observations, we make 50 update steps on all the data. As a result, we keep on updating the parameters as we explore more. We use Adam Optimizer with learning rate 10 −3 . At each step, we update over observations from all the arms. The flow of the exploration and parameter update is also shown in Algorithm 1.
We use the feature representations of the network used to model joint potentials described in Section 3.5 to also initialize the deep kernel induced by V. A final new linear layer of default output size 20 is added to project the feature representations.
In our proposed method BetaGP-SLP, described in Section 3.3, we take the kernel average of three neighbours for any arm with fewer than five observations.

D Simple Setting
In Section 3.2, we describe how the objective of BetaGP does not supervise the scale parameter. Further, in Section 3.3, we posit that the presence of sparse observations leads to learning a nonsmooth kernel. In this section, we illustrate these two observations using a simple setting.
We consider a simple estimation problem with 10 arms, their true accuracies go from 0.1 to a large value of 0.9 and then back to a small value of 0.2 as shown in the Table 8. In Table 8, we also show the number of observations per arm; Observe that the first three and the last three arms are sparsely observed.
We now present the fitted values by some of the methods we discussed in the main section. The index of an arm is the input for any estimator with no feature learning. Our motivation for discussing the simple setting is to illustrate the two limitations we discussed in the main content regarding the BetaGP objective: (a) the scale parameter of the BetaGP objective is not supervised (b) sparse observations lead to non-smooth surface. Toward these ends, we evaluate BetaGP, BetaGP-SL, BetaGP-SLP  Table 8: Arms, their indices, accuracies and number of observations (N) in the simple setting are shown in first three columns in that order. The scale parameter estimated using one of the algorithms for each arm is shown in the last three columns. Notice that BetaGP fitted scale parameter does not reflect the underlying observation sparsity for the first and last three arms.
methods on this setting. The fitted scale parameters for each arm by each of the estimators is shown in Table 8. Observe that BetaGP fitted scale parameter does not reflect the underlying observation sparsity of the first, last three arms. BetaGP-SL, BetaGP-SLP fitted scale values more faithfully reflect the underlying number of observations. All the numbers reported here are averaged over 20 seed runs.  Table 9: Bias-variance decomposition of MSE in the simple setting

Method
In Table 9, we show the bias 2 , variance decomposition of the mean squared error from 20 independent runs. Observe that BetaGP, BetaGP-SL have low bias but large variance and BetaGP-SLP has much lower variance at a slight expense of bias, as a result the overall MSE value for BetaGP-SLP is much lower than the other two. Moreover, we look at the fitted kernel length parameter (recall from Equation (7)) as a proxy for smoothness of the fitted kernel. Large kernel length is indicative of long range smoothness. The average kernel length for BetaGP, BetaGP-SL, BetaGP-SLP are 0.67, 0.68, 1.87 respectively. Despite using a GP kernel, we find the estimates of BetaGP, BetaGP-SL of high variance, that is also indicative of short range smoothness apparent from the low average kernel length. On the other hand, BetaGP-SLP imposes long range smoothness, as a result, decreases the MSE value more effectively when compared with BetaGP-SL.

E Calibration Training Details
In this section, we give further training details on the noise calibration method discussed in Section 3.5.
As discussed in Section 3.5, we use both labeled, small D and large U for training calibration parameters that are expressed in the objective (14). On the unlabeled data U , we use attribute values predicted using the predictors: {M k | k ∈ A} as a proxy for true values. The use of predicted values as the replacement for true value under-represents the attribute prediction error rate and interferes in the estimation of temperature parameters t. However, if we use U for training, we see a lot more attribute combinations and this can help identify more natural attribute combinations aiding in the learning of joint potential parameters of N . We mitigate the temperature estimation problem by up-sampling instances in D such that the loss in every batch contains equal contribution from D and U .
Recall that the MLE objective (14), contains contribution from two terms: (a) temperature scaled logits (b) attribute combination potential. In practice, we found that the second term (b) dominates the first, this causes under-training of the temperature parameters. Ideally, the two terms should be comparable and replaceable. We address this issue by dropping the second term corresponding to the network-assigned edge potential term in the objective half the times, which estimates better the temperature parameters. Further, we use a small held out fraction of D for network architecture search on L, and for early stopping. The training procedure is summarized in Alg. 2.

F More Task and Dataset Details
F.1 MF-CelebA, MF-IMDB As mentioned in Section 4.1, we hand-picked 12 binary attributes relevant for gender classification of the 40 total available attributes in the CelebA dataset. The twelve binary attributes are listed in Table 10, these constitute the A. A is the combination of twelve binary attributes and is 2 12 = 4, 096 large. The attributes related to hair color are retained in this list due to the recent finding that hair-color is spuriously correlated with the gender in CelebA [3]. We ignored several other gender-neutral or rare attributes. Chubby? 2 7 Mustache? 2 8 No Beard? 2 9 Wearing Hat? 2 10 Blurry? 2 11 Young? 2 12 Eyeglasses? 2 Table 10: Attribute list of MF-CelebA, MF-IMDB.

F.2 AC-COCOS, AC-COCOS10K
COCOS is a scene classification dataset, where pixel level supervision is provided. Methods are usually evaluated on pixel level classification accuracy. For simplicity, we cast it in to an object recognition task. The subset of ten animal labels we consider is shown in Table 11. We consider in our task five coarse background (stuff) labels by collapsing the fine labels to coarse using the label hierarchy shown in Table 12 and is the same as the official hierarchy [18].
We now describe how we cast the scene label classification task to an animal classification task. We first identify the subset of images in the train and validation set of the COCOS dataset which contain only one animal label, it could contain multiple background labels. If the image contains multiple animals, we exclude it, leaving around 23,000 images in the dataset. This is implemented in the routine: filter_ids_with_single_object of cocos3.py in the attached code. We also retrofit the scene classification models for animal classification. When the model (service model) labels pixels with more than one animal label, we retain the label associated with the largest number of pixels. This is implemented in fetch_preds routine of cocos3.py. Recall our calibration method makes use of the logit scores given by the attribute predictors, since we are aggregating prediction from multiple pixels, we do not have access to the logit scores. We simply set the logit score to +1 if a label is found in the image and -1 otherwise.
We follow the same protocol for both the tasks: AC-COCOS, AC-COCOS10K. The only difference between the two is the service model, AC-COCOS is a stronger service model trained on 164K size training data compared to AC-COCOS10K which is a model trained on a previous version of the dataset that is only 10K large. In these tasks, we use the same model for predicting the attributes and task labels since the pre-trained model we use is a scene label classifier. Both the pretrained 7 models were trained using ResNet101 architecture. Name bird cat dog horse sheep cow elephant bear zebra giraffe

G Statistics of Accuracy Surface
We show in Table 13, details about our two data sets such as the number of attributes, number of arms and number of active arms. Active arms are arms with a support of at least five and are the ones used for evaluation. Large number of arms as shown in the table exclude the possibility of manual supervision, since it is hard to obtain and label data that covers all the arms.
In Figures 14a and 14b, we show some statistics that illustrate the shape of the accuracy surface. We note that, although the service model's mean accuracy is high, the accuracy of the arms in the 10% quantile is abysmally low while arms in the top-quantiles have near perfect accuracy. This further motivates for why we need an accuracy surface instead of single accuracy estimate.
Dataset # attributes # arms # active arms CelebA 12 4096 398 COCOS 6 320 176 Table 13: Attribute statistics per dataset. First and second column show number of attributes and total possible combinations of the attributes. Third column shows number of attribute combinations (arms) with at least a support of five in the unlabeled data. These are the arms on which accuracy surface is evaluated.

H Standard Deviation and More Evaluation Metrics
We include standard deviation accompanying numbers in Table 2 for macro MSE in Table 15 and for worst MSE in Table 17.
In the main content of the paper, we gave results using macro and worst MSE. In this following section, we show results using two other metrics. We follow the same setup as in Section 4.5.

Micro-averaged MSE:
We assign importance to each arm based on its support. The error per arm is multiplied by its support (in U ). Results with this error are shown in Table 16. The best predictor with this metric is the point estimate given by the CPredictor estimator which is not surprising since very few arms with high frequency dominate this metric. Infrequent MSE: In Table 18, we show MSE evaluated only on the 50 arms that are least frequent in U . For each of the above evaluation metrics, the trend between BetaGP BetaGP-SL, BetaGP-SLP is statistically significant.

I Pseudocode
The full flow of estimation and exploration is summarized in Algorithm 1. Algorithm 2 shows the calibration training sub-routine.   Table 16: Micro-averaged MSE along with standard deviation on all tasks. Shown after trailing '/' is the standard deviation.