summaries: This paper considers the robust estimation problem under Huber’s \epsilon-contamination model. This problem is a hot topic in theoretical statistics and theoretical computer science community in recent 3 years.  From theoretical statistics community, the main approach is through depth functions. Solving the robust estimation problem can be reduced to solving a min-max problem. While the formulation is clean and can achieve the optimal statistical rate, solving the min-max problem is computationally intractable in general. On the other hand, approaches from TCS community are more involved and sometimes cannot achieve the optimal statistical rate (especially for the general distribution).

summaries: The paper considers the problem of robust high dimensional estimation in Huber’s contamination model. The algorithm is given samples from a distribution (1 - eps) * P + eps * Q, where P is a “nice” distribution (e.g. a Gaussian), eps is the fraction of contaminated points, and Q is some unconstrained noise distribution. The goal is then to estimate parameters of P as well as possible, given this noise. The settings they primarily consider in this paper are when P is a Gaussian with unknown mean and identity covariance, or when it is a Gaussian with unknown covariance. Classical estimators such as Tukey depth or matrix depth for these problems achieve optimal minimax rates, but are computationally expensive to compute. However, recent work of [1,2] propose efficient estimators for this problem that (nearly) achieve these rates.

abstractText: Robust estimation under Huber’s -contamination model has become an important topic in statistics and theoretical computer science. Statistically optimal procedures such as Tukey’s median and other estimators based on depth functions are impractical because of their computational intractability. In this paper, we establish an intriguing connection between f -GANs and various depth functions through the lens of f -Learning. Similar to the derivation of f GANs, we show that these depth functions that lead to statistically optimal robust estimators can all be viewed as variational lower bounds of the total variation distance in the framework of f -Learning. This connection opens the door of computing robust estimators using tools developed for training GANs. In particular, we show in both theory and experiments that some appropriate structures of discriminator networks with hidden layers in GANs lead to statistically optimal robust location estimators for both Gaussian distribution and general elliptical distributions where first moment may not exist.

authors:
Syed Mumtaz Ali
Samuel D Silvey

title: A general class of coefficients of divergence of one distribution from another

venue: Journal of the Royal Statistical Society. Series B (Methodological),

authors:
Nina Amenta
Marshall Bern
David Eppstein
S-H Teng

title: Regression depth and center points

venue: Discrete & Computational Geometry,

title: For valid generalization the size of the weights is more important than the size of the network

venue: In Advances in neural information processing systems,

authors:
Peter L Bartlett
Shahar Mendelson

title: Rademacher and gaussian complexities: Risk bounds and structural results

venue: Journal of Machine Learning Research,

title: An optimal randomized algorithm for maximum tukey depth

venue: In Proceedings of the fifteenth annual ACM-SIAM symposium on Discrete algorithms,

authors:
Mengjie Chen
Chao Gao
Zhao Ren

title: Robust covariance and scatter matrix estimation under hubers contamination model

venue: The Annals of Statistics,

title: Eine informationstheoretische ungleichung und ihre anwendung auf beweis der ergodizitaet von markoffschen ketten

venue: Magyer Tud. Akad. Mat. Kutato Int. Koezl.,

authors:
Luc Devroye
Gábor Lugosi

title: Combinatorial methods in density estimation

authors:
Ilias Diakonikolas
Gautam Kamath
Daniel M Kane
Jerry Li
Ankur Moitra
Alistair Stewart

title: Robust estimators in high dimensions without the computational intractability

venue: In Foundations of Computer Science (FOCS),

authors:
Ilias Diakonikolas
Daniel Kane
Alistair Stewart

title: Robust learning of fixed-structure bayesian networks

venue: arXiv preprint arXiv:1606.07384,

title: Being robust (in high dimensions) can be practical

venue: arXiv preprint arXiv:1703.00893,

authors:
Ilias Diakonikolas
Gautam Kamath
Daniel M Kane
Jerry Li
Jacob Steinhardt
Alistair Stewart

title: Sever: A robust meta-algorithm for stochastic optimization

venue: arXiv preprint arXiv:1803.02815,

authors:
Ilias Diakonikolas
Daniel M Kane
Alistair Stewart

title: List-decodable robust mean estimation and learning mixtures of spherical gaussians

venue: In Proceedings of the 50th Annual ACM SIGACT Symposium on Theory of Computing,

authors:
Ilias Diakonikolas
Weihao Kong
Alistair Stewart

title: Efficient algorithms and lower bounds for robust linear regression

venue: arXiv preprint arXiv:1806.00040,

authors:
David L Donoho
Miriam Gasko

title: Breakdown properties of location estimates based on halfspace depth and projected outlyingness

authors:
Simon S Du
Sivaraman Balakrishnan
Aarti Singh

title: Computationally efficient robust estimation of sparse functionals

venue: arXiv preprint arXiv:1702.07709,

title: Symmetric Multivariate and Related Distributions: 0

title: Robust regression via mutivariate regression depth

venue: arXiv preprint arXiv:1702.04656,

authors:
Xavier Glorot
Yoshua Bengio

title: Understanding the difficulty of training deep feedforward neural networks

venue: In Proceedings of the thirteenth international conference on artificial intelligence and statistics,

authors:
Ian Goodfellow
Jean Pouget-Abadie
Mehdi Mirza
Bing Xu
David Warde-Farley
Sherjil Ozair
Aaron Courville
Yoshua Bengio

title: Generative adversarial nets

title: Robust estimation of a location parameter

venue: The annals of mathematical statistics,

title: A robust version of the probability ratio test

venue: The Annals of Mathematical Statistics,

authors:
Pravesh K Kothari
Jacob Steinhardt
David Steurer

title: Robust moment estimation and improved clustering via sum of squares

authors:
Karol Kurach
Mario Lucic
Xiaohua Zhai
Marcin Michalski
Sylvain Gelly

title: The gan landscape: Losses, architectures, regularization, and normalization

venue: arXiv preprint arXiv:1807.04720,

authors:
Kevin A Lai
Anup B Rao
Santosh Vempala

title: Agnostic estimation of mean and covariance

authors:
Shuang Liu
Olivier Bousquet
Kamalika Chaudhuri

title: Approximation and convergence properties of generative adversarial learning

authors:
Mario Lucic
Karol Kurach
Marcin Michalski
Sylvain Gelly
Olivier Bousquet

title: Are gans created equal? a large-scale study

venue: arXiv preprint arXiv:1711.10337,

title: On the method of bounded differences

venue: Surveys in combinatorics,

authors:
Takeru Miyato
Toshiki Kataoka
Masanori Koyama
Yuichi Yoshida

title: Spectral normalization for generative adversarial networks

venue: arXiv preprint arXiv:1802.05957,

title: On depth and deep points: a calculus

authors:
Ivan Mizera
Christine H Müller

title: Location–scale depth

venue: Journal of the American Statistical Association,

authors:
Youssef Mroueh
Tom Sercu
Vaibhava Goel

title: Mcgan: Mean and covariance feature matching gan

venue: arXiv preprint arXiv:1702.08398,

authors:
XuanLong Nguyen
Martin J Wainwright
Michael I Jordan

title: Estimating divergence functionals and the likelihood ratio by convex risk minimization

venue: IEEE Transactions on Information Theory,

authors:
Sebastian Nowozin
Botond Cseke
Ryota Tomioka

title: f-gan: Training generative neural samplers using variational divergence minimization

authors:
Davy Paindaveine
Germain Van Bever

title: Halfspace depths for scatter, concentration and shape matrices

venue: arXiv preprint arXiv:1704.06160,

title: Convergence of stochastic processes

authors:
Alec Radford
Luke Metz
Soumith Chintala

title: Unsupervised representation learning with deep convolutional generative adversarial networks

venue: arXiv preprint arXiv:1511.06434,

authors:
Peter J Rousseeuw
Mia Hubert

title: Regression depth

authors:
Peter J Rousseeuw
Anja Struyf

title: Computing location depth and regression depth in higher dimensions

venue: Statistics and Computing,

authors:
Tim Salimans
Ian Goodfellow
Wojciech Zaremba
Vicki Cheung
Alec Radford
Xi Chen

title: Improved techniques for training gans

title: Mathematics and the picturing of data

venue: In Proceedings of the International Congress of Mathematicians, Vancouver,

authors:
Marc van Kreveld
Joseph SB Mitchell
Peter Rousseeuw
Micha Sharir
Jack Snoeyink
Bettina Speckmann

title: Efficient algorithms for maximum regression depth

venue: In Proceedings of the fifteenth annual symposium on Computational geometry,

title: Rates of convergence of minimum distance estimators and kolmogorov’s entropy

title: Some extensions of tukey’s depth function

venue: Journal of Multivariate Analysis,

title: 2018) is also worth trying since it can prevent the weight from going to infinity and

heading: 1 INTRODUCTION

text:
In the setting of Huber’s -contamination model (Huber, 1964; 1965), one has i.i.d observations
X1, ..., Xn ∼ (1− )Pθ + Q, (1) and the goal is to estimate the model parameter θ. Under the data generating process (1), each observation has a 1 − probability to be drawn from Pθ and the other probability to be drawn from the contamination distributionQ. The presence of an unknown contamination distribution poses both statistical and computational challenges to the problem. For example, consider a normal mean estimation problem with Pθ = N(θ, Ip). Due to the contamination of data, the sample average, which is optimal when = 0, can be arbitrarily far away from the true mean if Q charges a positive probability at infinity. Moreover, even robust estimators such as coordinatewise median and geometric median are proved to be suboptimal under the setting of (1) (Chen et al., 2018; Diakonikolas et al., 2016a; Lai et al., 2016). The search for both statistically optimal and computationally feasible procedures has become a fundamental problem in areas including statistics and computer science.
For the normal mean estimation problem, it has been shown in Chen et al. (2018) that the minimax rate with respect to the squared `2 loss is pn ∨
2, and is achieved by Tukey’s median (Tukey, 1975). Despite the statistical optimality of Tukey’s median, its computation is not tractable. In fact, even an approximate algorithm takes O(eCp) in time (Amenta et al., 2000; Chan, 2004; Rousseeuw & Struyf, 1998).
Recent developments in theoretical computer science are focused on the search of computationally tractable algorithms for estimating θ under Huber’s -contamination model (1). The success of the efforts started from two fundamental papers Diakonikolas et al. (2016a); Lai et al. (2016), where two different but related computational strategies “iterative filtering” and “dimension halving” were proposed to robustly estimate the normal mean. These algorithms can provably achieve the minimax rate pn ∨
2 up to a poly-logarithmic factor in polynomial time. The main idea behind the two methods is a critical fact that a good robust moment estimator can be certified efficiently by higher moments. This idea was later further extended (Diakonikolas et al., 2017; Du et al., 2017; Diakonikolas et al., 2016b; 2018a;c;b; Kothari et al., 2018) to develop robust and computable procedures for various other problems.
However, many of the computationally feasible procedures for robust mean estimation in the literature rely on the knowledge of covariance matrix and sometimes the knowledge of contamination proportion. Even though these assumptions can be relaxed, nontrivial modifications of the algorithms are required for such extensions and statistical error rates may also be affected. Compared with these computationally feasible procedures proposed in the recent literature for robust estimation, Tukey’s median (9) and other depth-based estimators (Rousseeuw & Hubert, 1999; Mizera, 2002; Zhang, 2002; Mizera & Müller, 2004; Paindaveine & Van Bever, 2017) have some indispensable advantages in terms of their statistical properties. First, the depth-based estimators have clear objective functions that can be interpreted from the perspective of projection pursuit (Mizera, 2002). Second, the depth-based procedures are adaptive to unknown nuisance parameters in the models such as covariance structures, contamination proportion, and error distributions (Chen et al., 2018; Gao, 2017). Last but not least, Tukey’s depth and other depth functions are mostly designed for robust quantile estimation, while the recent advancements in the theoretical computer science literature are all focused on robust moments estimation. Although this is not an issue when it comes to normal mean estimation, the difference is fundamental for robust estimation under general settings such as elliptical distributions where moments do not necessarily exist.
Given the desirable statistical properties discussed above, this paper is focused on the development of computational strategies of depth-like procedures. Our key observation is that robust estimators that are maximizers of depth functions, including halfspace depth, regression depth and covariance matrix depth, can all be derived under the framework of f -GAN (Nowozin et al., 2016). As a result, these depth-based estimators can be viewed as minimizers of variational lower bounds of the total variation distance between the empirical measure and the model distribution (Proposition 2.1). This observation allows us to leverage the recent developments in the deep learning literature to compute these variational lower bounds through neural network approximations. Our theoretical results give insights on how to choose appropriate neural network classes that lead to minimax optimal robust estimation under Huber’s -contamination model. In particular, Theorem 3.1 and 3.2 characterize the networks which can robustly estimate the Gaussian mean by TV-GAN and JS-GAN, respectively; Theorem 4.1 is an extension to robust location estimation under the class of elliptical distributions which includes Cauchy distribution whose mean does not exist. Numerical experiments in Section 5 are provided to show the success of these GANs.
2 ROBUST ESTIMATION AND f -GAN We start with the definition of f -divergence (Csiszár, 1964; Ali & Silvey, 1966). Given a strictly convex function f that satisfies f(1) = 0, the f -GAN between two probability distributions P and Q is defined by
Df (P‖Q) = ∫ f ( p
q
) dQ. (2)
Here, we use p(·) and q(·) to stand for the density functions of P and Q with respect to some common dominating measure. For a fully rigorous definition, see Polyanskiy & Wu (2017). Let f∗ be the convex conjugate of f . That is, f∗(t) = supu∈domf (ut− f(u)). A variational lower bound of (2) is
Df (P‖Q) ≥ sup T∈T [EPT (X)− EQf∗(T (X))] . (3)
Note that the inequality (3) holds for any class T , and it becomes an equality whenever the class T contains the function f ′ (p/q) (Nguyen et al., 2010). For notational simplicity, we also use f ′ for an arbitrary element of
the subdifferential when the derivative does not exist. With i.i.d. observations X1, ..., Xn ∼ P , the variational lower bound (3) naturally leads to the following learning method
P̂ = argmin Q∈Q sup T∈T
[ 1
n n∑ i=1 T (Xi)− EQf∗(T (X))
] . (4)
The formula (4) is a powerful and general way to learn the distribution P from its i.i.d. observations. It is known as f -GAN (Nowozin et al., 2016), an extension of GAN (Goodfellow et al., 2014), which stands for generative adversarial networks. The idea is to find a P̂ so that the best discriminator T in the class T cannot tell the difference between P̂ and the empirical distribution 1n ∑n i=1 δXi .
2.1 f -LEARNING: A UNIFIED FRAMEWORK
Our f -Learning framework is based on a special case of the variational lower bound (3). That is,
Df (P‖Q) ≥ sup Q̃∈Q̃Q
[ EP f ′ ( q̃(X)
q(X)
) − EQf∗ ( f ′ ( q̃(X)
q(X)
))] , (5)
where q̃(·) stands for the density function of Q̃. Note that here we allow the class Q̃Q to depend on the distribution Q in the second argument of Df (P‖Q). Compare (5) with (3), and it is easy to realize that (5) is a special case of (3) with
T = TQ = { f ′ ( q̃
q
) : q̃ ∈ Q̃Q } . (6)
Moreover, the inequality (5) becomes an equality as long as P ∈ Q̃Q. The sample version of (5) leads to the following learning method
P̂ = argmin Q∈Q sup Q̃∈Q̃Q
[ 1
n n∑ i=1 f ′ ( q̃(Xi) q(Xi) ) − EQf∗ ( f ′ ( q̃(X) q(X) ))] . (7)
The learning method (7) will be referred to as f -Learning in the sequel. It is a very general framework that covers many important learning procedures as special cases. For example, consider the special case where Q̃Q = Q̃ independent of Q, Q = Q̃, and f(x) = x log x. Direct calculations give f ′(x) = log x + 1 and f∗(t) = et−1. Therefore, (7) becomes
P̂ = argmin Q∈Q sup Q̃∈Q
1
n n∑ i=1 log q̃(Xi) q(Xi) = argmax q∈Q 1 n n∑ i=1 log q(Xi),
which is the maximum likelihood estimator (MLE).

heading: 2.2 TV-LEARNING AND DEPTH-BASED ESTIMATORS

text:
An important generator f that we will discuss here is f(x) = (x−1)+. This leads to the total variation distance Df (P‖Q) = 12 ∫ |p− q|. With f ′(x) = I{x ≥ 1} and f∗(t) = tI{0 ≤ t ≤ 1}, the TV-Learning is given by
P̂ = argmin Q∈Q sup Q̃∈QQ
[ 1
n n∑ i=1 I { q̃(Xi) q(Xi) ≥ 1 } −Q ( q̃ q ≥ 1 )] . (8)
A closely related idea was previously explored by Yatracos (1985); Devroye & Lugosi (2012). The following proposition shows that when Q̃Q approaches toQ in some neighborhood, TV-Learning leads to robust estimators that are defined as the maximizers of various depth functions including Tukey’s depth, regression depth, and covariance depth.
Proposition 2.1. The TV-Learning (8) includes the following special cases:
1. Tukey’s halfspace depth: Take Q = {N(η, Ip) : η ∈ Rp} and Q̃η = {N(η̃, Ip) : ‖η̃ − η‖ ≤ r}. As r → 0, (8) becomes
θ̂ = argmax η∈Rp inf ‖u‖=1
1
n n∑ i=1 I{uT (Xi − η) ≥ 0}. (9)
2. Regression depth: Take Q = { Py,X = Py|XPX : Py|X = N(X T η, 1), η ∈ Rp }
, and Q̃η = { Py,X = Py|XPX : Py|X = N(X T η̃, 1), ‖η̃ − η‖ ≤ r } . As r → 0, (8) becomes
θ̂ = argmax η∈Rp inf ‖u‖=1
1
n n∑ i=1 I{uTXi(yi −XTi η) ≥ 0}. (10)
3. Covariance matrix depth: Take Q = {N(0,Γ) : Γ ∈ Ep}, where Ep stands for the class of p × p covariance matrices, and Q̃Γ = { N(0, Γ̃) : Γ̃−1 = Γ−1 + r̃uuT ∈ Ep, |r̃| ≤ r, ‖u‖ = 1 } . As r → 0,
(8) becomes
Σ̂ = argmin Γ∈Ep sup ‖u‖=1
[( 1
n n∑ i=1 I{|uTXi|2 ≤ uTΓu} − P(χ21 ≤ 1)
) (11)
∨
( 1
n n∑ i=1 I{|uTXi|2 > uTΓu} − P(χ21 > 1)
)] .
The formula (9) is recognized as Tukey’s median, the maximizer of Tukey’s halfspace depth. A traditional understanding of Tukey’s median is that (9) maximizes the halfspace depth (Donoho & Gasko, 1992) so that θ̂ is close to the centers of all one-dimensional projections of the data. In the f -Learning framework, N(θ̂, Ip) is understood to be the minimizer of a variational lower bound of the total variation distance. The formula (10) gives the estimator that maximizes the regression depth proposed by Rousseeuw & Hubert (1999). It is worth noting that the derivation of (10) does not depend on the marginal distribution PX in the linear regression model. Finally, (11) is related to the covariance matrix depth (Zhang, 2002; Chen et al., 2018; Paindaveine & Van Bever, 2017). All of the estimators (9), (10) and (11) are proved to achieve the minimax rate for the corresponding problems under Huber’s -contamination model (Chen et al., 2018; Gao, 2017).
2.3 FROM f -LEARNING TO f -GAN
The connection to various depth functions shows the importance of TV-Learning in robust estimation. However, it is well-known that depth-based estimators are very hard to compute (Amenta et al., 2000; van Kreveld et al., 1999; Rousseeuw & Struyf, 1998), which limits their applications only for very low-dimensional problems. On the other hand, the general f -GAN framework (4) has been successfully applied to learn complex distributions and images in practice (Goodfellow et al., 2014; Radford et al., 2015; Salimans et al., 2016). The major difference that gives the computational advantage to f -GAN is its flexibility in terms of designing the discriminator class T using neural networks compared with the pre-specified choice (6) in f -Learning. While f -Learning provides a unified perspective in understanding various depth-based procedures in robust estimation, we can step back into the more general f -GAN for its computational advantages, and to design efficient computational strategies.

heading: 3 ROBUST MEAN ESTIMATION VIA GAN

text: In this section, we focus on the problem of robust mean estimation under Huber’s -contamination model. Our goal is to reveal how the choice of the class of discriminators affects robustness and statistical optimality under the simplest possible setting. That is, we have i.i.d. observations X1, ..., Xn ∼ (1− )N(θ, Ip) + Q, and we need to estimate the unknown location θ ∈ Rp with the contaminated data. Our goal is to achieve the minimax rate pn ∨ 2 with respect to the squared `2 loss uniformly over all θ ∈ Rp and all Q.

heading: 3.1 RESULTS FOR TV-GAN

text:
We start with the total variation GAN (TV-GAN) with f(x) = (x − 1)+ in (4). For the Gaussian location family, (4) can be written as
θ̂ = argmin η∈Rp max D∈D
[ 1
n n∑ i=1 D(Xi)− EN(η,Ip)D(X)
] , (12)
with T (x) = D(x) in (4). Now we need to specify the class of discriminatorsD to solve the classification problem between N(η, Ip) and the empirical distribution 1n ∑n i=1 δXi . One of the simplest discriminator classes is the logistic regression, D = { D(x) = sigmoid(wTx+ b) : w ∈ Rp, b ∈ R } . (13)
With D(x) = sigmoid(wTx+ b) = (1 + e−w T x−b)−1 in (13), the procedure (12) can be viewed as a smoothed version of TV-Learning (8). To be specific, the sigmoid function sigmoid(wTx + b) tends to an indicator function as ‖w‖ → ∞, which leads to a procedure very similar to (9). In fact, the class (13) is richer than the one used in (9), and thus (12) can be understood as the minimizer of a sharper variational lower bound than that of (9).
Theorem 3.1. Assume pn + 2 ≤ c for some sufficiently small constant c > 0. With i.i.d. observations X1, ..., Xn ∼ (1− )N(θ, Ip) + Q, the estimator θ̂ defined by (12) satisfies
‖θ̂ − θ‖2 ≤ C ( p n ∨ 2 ) ,
with probability at least 1 − e−C′(p+n 2) uniformly over all θ ∈ Rp and all Q. The constants C,C ′ > 0 are universal.
Though TV-GAN can achieve the minimax rate pn ∨ 2 under Huber’s contamination model, it may suffer from optimization difficulties especially when the distributions Q and N(θ, Ip) are far away from each other, as shown in Figure 1.

heading: 3.2 RESULTS FOR JS-GAN

text:
Given the intractable optimization property of TV-GAN, we next turn to Jensen-Shannon GAN (JS-GAN) with f(x) = x log x− (x+ 1) log x+12 . The estimator is defined by
θ̂ = argmin η∈Rp max D∈D
[ 1
n n∑ i=1 logD(Xi) + EN(η,Ip) log(1−D(Xi))
] + log 4, (14)
with T (x) = logD(x) in (4). This is exactly the original GAN (Goodfellow et al., 2014) specialized to the normal mean estimation problem. The advantages of JS-GAN over other forms of GAN have been studied extensively in the literature (Lucic et al., 2017; Kurach et al., 2018).
Unlike TV-GAN, our experiment results show that (14) with the logistic regression discriminator class (13) is not robust to contamination. However, if we replace (13) by a neural network class with one or more hidden layers, the estimator will be robust and will also work very well numerically.
To understand why and how the class of the discriminators affects the robustness property of JS-GAN, we introduce a new concept called restricted Jensen-Shannon divergence. Let g : Rp → Rd be a function that maps a p-dimensional observation to a d-dimensional feature space. The restricted Jensen-Shannon divergence between two probability distributions P and Q with respect to the feature g is defined as
JSg(P,Q) = max w∈W
[ EP log sigmoid(w T g(X)) + EQ log(1− sigmoid(wT g(X))) ] + log 4.
In other words, P and Q are distinguished by a logistic regression classifier that uses the feature g(X). It is easy to see that JSg(P,Q) is a variational lower bound of the original Jensen-Shannon divergence. The key property of JSg(P,Q) is given by the following proposition.
Proposition 3.1. AssumeW is a convex set that contains an open neighborhood of 0. Then, JSg(P,Q) = 0 if and only if EP g(X) = EQg(X).
The proposition asserts that JSg(·, ·) cannot distinguish P and Q if the feature g(X) has the same expected value under the two distributions. This generalized moment matching effect has also been studied by Liu et al. (2017) for general f -GANs. However, the linear discriminator class considered in Liu et al. (2017) is parameterized in a different way compared with the discriminator class here.
When we apply Proposition 3.1 to robust mean estimation, the JS-GAN is trying to match the values of 1 n ∑n i=1 g(Xi) and EN(η,Iη)g(X) for the feature g(X) used in the logistic regression classifier. This explains what we observed in our numerical experiments. A neural net without any hidden layer is equivalent to a logistic regression with a linear feature g(X) = (XT , 1)T ∈ Rp+1. Therefore, whenever η = 1n ∑n i=1Xi, we
have JSg ( 1 n ∑n i=1 δXi , N(η, Ip) ) = 0, which implies that the sample mean is a global maximizer of (14). On the other hand, a neural net with at least one hidden layers involves a nonlinear feature function g(X), which is the key that leads to the robustness of (14).
We will show rigorously that a neural net with one hidden layer is sufficient to make (14) robust and optimal. Consider the following class of discriminators,
D = D(x) = sigmoid ∑ j≥1 wjσ(u T j x+ bj)  : ∑ j≥1 |wj | ≤ κ, uj ∈ Rp, bj ∈ R  . (15) The class (15) consists of two-layer neural network functions. While the dimension of the input layer is p, the dimension of the hidden layer can be arbitrary, as long as the weights have a bounded `1 norm. The nonlinear activation function σ(·) is allowed to take 1) indicator: σ(x) = I{x ≥ 1}, 2) sigmoid: σ(x) = 11+e−x , 3) ramp: σ(x) = max(min(x + 1/2, 1), 0). Other bounded activation functions are also possible, but we do not exclusively list them. The rectified linear unit (ReLU) will be studied in Appendix A.
Theorem 3.2. Consider the estimator θ̂ defined by (14) with D specified by (15). Assume pn + 2 ≤ c for some sufficiently small constant c > 0, and set κ = O (√
p n + ) . With i.i.d. observations X1, ..., Xn ∼
(1− )N(θ, Ip) + Q, we have ‖θ̂ − θ‖2 ≤ C ( p n ∨ 2 ) , with probability at least 1 − e−C′(p+n 2) uniformly over all θ ∈ Rp and all Q. The constants C,C ′ > 0 are universal.

heading: 4 ELLIPTICAL DISTRIBUTIONS

text:
An advantage of Tukey’s median (9) is that it leads to optimal robust location estimation under general elliptical distributions such as Cauchy distribution whose mean does not exist. In this section, we show that JS-GAN shares the same property. A random vectorX ∈ Rp follows an elliptical distribution if it admits a representation
X = θ + ξAU,
where U is uniformly distributed on the unit sphere {u ∈ Rp : ‖u‖ = 1} and ξ ≥ 0 is a random variable independent of U that determines the shape of the elliptical distribution (Fang, 2017). The center and the scatter matrix is θ and Σ = AAT .
For a unit vector v, let the density function of ξvTU be h. Note that h is independent of v because of the symmetry of U . Then, there is a one-to-one relation between the distribution of ξ and h, and thus the triplet (θ,Σ, h) fully parametrizes an elliptical distribution.
Note that h and Σ = AAT are not identifiable, because ξA = (cξ)(c−1A) for any c > 0. Therefore, without loss of generality, we can restrict h to be a member of the following class
H = { h : h(t) = h(−t), h ≥ 0, ∫ h = 1, ∫ σ(t)(1− σ(t))h(t)dt = 1 } .
This makes the parametrization (θ,Σ, h) of an elliptical distribution fully identifiable, and we use EC(θ,Σ, h) to denote an elliptical distribution parametrized in this way.
The JS-GAN estimator is defined as
(θ̂, Σ̂, ĥ) = argmin η∈Rp,Γ∈Ep(M),g∈H max D∈D
[ 1
n n∑ i=1 logD(Xi) + EEC(η,Γ,g) log(1−D(X))
] + log 4, (16)
where Ep(M) is the set of all positive semi-definite matrix with spectral norm bounded by M .
Theorem 4.1. Consider the estimator θ̂ defined above withD specified by (15). AssumeM = O(1), pn+ 2 ≤ c for some sufficiently small constant c > 0, and set κ = O (√
p n + ) . With i.i.d. observations X1, ..., Xn ∼
(1− )EC(θ,Σ, h) + Q, we have ‖θ̂ − θ‖2 ≤ C ( p n ∨ 2 ) , with probability at least 1 − e−C′(p+n 2) uniformly over all θ ∈ Rp, Σ ∈ Ep(M) and all Q. The constants C,C ′ > 0 are universal.
Remark 4.1. The result of Theorem 4.1 also holds (and is proved) under the strong contamination model (Diakonikolas et al., 2016a). That is, we have i.i.d. observations X1, ..., Xn ∼ P for some P satisfying TV(P,EC(θ,Σ, h)) ≤ . See its proof in Appendix D.2.
Note that Theorem 4.1 guarantees the same convergence rate as in the Gaussian case for all elliptical distributions. This even includes multivariate Cauchy where mean does not exist. Therefore, the location estimator (16) is fundamentally different from Diakonikolas et al. (2016a); Lai et al. (2016), which is only designed for robust mean estimation. We will show such a difference in our numerical results.
To achieve rate-optimality for robust location estimation under general elliptical distributions, the estimator (16) is different from (14) only in the generator class. They share the same discriminator class (15). This underlines
an important principle for designing GAN estimators: the overall statistical complexity of the estimator is only determined by the discriminator class.
The estimator (16) also outputs (Σ̂, ĥ), but we do not claim any theoretical property for (Σ̂, ĥ) in this paper. This will be systematically studied in a future project.

heading: 5 NUMERICAL EXPERIMENTS

text: In this section, we give extensive numerical studies of robust mean estimation via GAN. After introducing the implementation details in Section 5.1, we verify our theoretical results on minimax estimation with both TVGAN and JS-GAN in Section 5.2. Comparison with other methods on robust mean estimation in the literature is given in Section 5.3. The effects of various network structures are studied in Section 5.4. Adaptation to unknown covariance is studied in Section 5.5. In all these cases, we assume i.i.d. observations are drawn from (1− )N(0p, Ip) + Q with and Q to be specified. Finally, adaptation to elliptical distributions is studied in Section 5.6.

heading: 5.1 IMPLEMENTATIONS

text: We adopt the standard algorithmic framework of f -GANs (Nowozin et al., 2016) for the implementation of JSGAN and TV-GAN for robust mean estimation. In particular, the generator for mean estimation is Gη(Z) = Z + η with Z ∼ N(0p, Ip); the discriminator D is a multilayer perceptron (MLP), where each layer consisting of a linear map and a sigmoid activation function and the number of nodes will vary in different experiments to be specified below. Details related to algorithms, tuning, critical hyper-parameters, structures of discriminator networks and other training tricks for stabilization and acceleration are discussed in Appendix B.1. A PyTorch implementation is available at https://github.com/zhuwzh/Robust-GAN-Center.

heading: 5.2 NUMERICAL SUPPORTS FOR THE MINIMAX RATES

text:
We verify the minimax rates achieved by TV-GAN (Theorem 3.1) and JS-GAN (Theorem 3.2) via numerical experiments. Two main scenarios we consider here are √ p/n < and √ p/n > , where in both cases, various types of contamination distributions Q are considered. Specifically, the choice of contamination distributions Q includes N(µ ∗ 1p, Ip) with µ ranges in {0.2, 0.5, 1, 5}, N(0.5 ∗ 1p,Σ) and Cauchy(τ ∗ 1p). Details of the construction of the covariance matrix Σ is given in Appendix B.2. The distribution Cauchy(τ ∗ 1p) is obtained by combining p independent one-dimensional standard Cauchy with location parameter τj = 0.5.
The main experimental results are summarized in Figure 2, where the `2 error we present is the maximum error among all choices of Q, and detailed numerical results can be founded in Tables 7, 8 and 9 in Appendix. We separately explore the relation between the error and one of , √ p and 1/ √ n with the other two parameters
fixed. The study of the relation between the `2 error and is in the regime √ p/n < so that dominates
the minimax rate. The scenario √ p/n > is considered in the study of the effects of √ p and 1/ √ n. As is shown in Figure 2, the errors are approximately linear against the corresponding parameters in all cases, which empirically verifies the conclusions of Theorem 3.1 and Theorem 3.2.

heading: 5.3 COMPARISONS WITH OTHER METHODS

text:
We perform additional experiments to compare with other methods including dimension halving (Lai et al., 2016) and iterative filtering (Diakonikolas et al., 2017) under various settings. We emphasize that our method does not require any knowledge about the nuisance parameters such as the contamination proportion . Tuning GAN is only a matter of optimization and one can tune parameters based on the objective function only.
Table 1 shows the performances of JS-GAN, TV-GAN, dimension halving, and iterative filtering. The network structure, for both JS-GAN and TV-GAN, has one hidden layer with 20 hidden units when the sample size is 50,000 and 2 hidden units when sample size is 5,000. The critical hyper-parameters we apply is given in Appendix and it turns out that the choice of the hyper-parameter is robust against different models when the net structures are the same. To summarize, our method outperforms other algorithms in most cases. TV-GAN is good at cases when Q and N(0p, Ip) are non-separable but fails when Q is far away from N(0p, Ip) due to optimization issues discussed in Section 3.1 (Figure 1). On the other hand, JS-GAN stably achieves the lowest error in separable cases and also shows competitive performances for non-separable ones.

heading: 5.4 NETWORK STRUCTURES

text: We further study the performance of JS-GAN with various structures of neural networks. The main observation is tuning networks with one-hidden layer becomes tough as the dimension grows (e.g. p ≥ 200), while a deeper network can significantly refine the situation perhaps by improving the landscape. Some experiment results are given in Table 2. On the other hand, one-hidden layer performs not worse than deeper networks when dimension is not very large (e.g. p ≤ 100). More experiments are given in Appendix B.4. Additional theoretical results for deep neural nets are given in Appendix A.

heading: 5.5 ADAPTATION TO UNKNOWN COVARIANCE

text:
The robust mean estimator constructed through JS-GAN can be easily made adaptive to unknown covariance structure, which is a special case of (16). We define
(θ̂, Σ̂) = argmin η∈Rp,Γ∈Ep max D∈D
[ 1
n n∑ i=1 logD(Xi) + EN(η,Γ) log(1−D(Xi))
] + log 4,
The estimator θ̂, as a result, is rate-optimal even when the true covariance matrix is not necessarily identity and is unknown (see Theorem 4.1). Below, we demonstrate some numerical evidence of the optimality of θ̂ as well as the error of Σ̂ in Table 3.

heading: 5.6 ADAPTATION TO ELLIPTICAL DISTRIBUTIONS

text:
We consider the estimation of the location parameter θ in elliptical distribution EC(θ,Σ, h) by the JS-GAN defined in (16). In particular, we study the case with i.i.d. observationsX1, ..., Xn ∼ (1− )Cauchy(θ, Ip)+ Q. The density function of Cauchy(θ,Σ) is given by p(x; θ,Σ) ∝ |Σ|−1/2 ( 1 + (x− θ)TΣ−1(x− θ) )−(1+p)/2 .
Compared with Algorithm (1), the difference lies in the choice of the generator. We consider the generator G1(ξ, U) = gω(ξ)U + θ, where gω(ξ) is a non-negative neural network parametrized by ω and some random variable ξ. The random vector U is sampled from the uniform distribution on {u ∈ Rp : ‖u‖ = 1}. If the scatter matrix is unknown, we will use the generatorG2(ξ, U) = gω(ξ)AU+θ, withAAT modeling the scatter matrix.
Table 4 shows the comparison with other methods. Our method still works well under Cauchy distribution, while the performance of other methods that rely on moment conditions deteriorates in this setting.

heading: ACKNOWLEDGEMENT

text: The research of Chao Gao was supported in part by NSF grant DMS-1712957 and NSF Career Award DMS1847590. The research of Yuan Yao was supported in part by Hong Kong Research Grant Council (HKRGC) grant 16303817, National Basic Research Program of China (No. 2015CB85600), National Natural Science Foundation of China (No. 61370004, 11421110001), as well as awards from Tencent AI Lab, Si Family Foundation, Baidu Big Data Institute, and Microsoft Research-Asia.

heading: A ADDITIONAL THEORETICAL RESULTS

text:
In this section, we investigate the performance of discriminator classes of deep neural nets with the ReLU activation function. Since our goal is to learn a p-dimensional mean vector, a deep neural network discriminator without any regularization will certainly lead to overfitting. Therefore, it is crucial to design a network class with some appropriate regularizations. Inspired by the work of Bartlett (1997); Bartlett & Mendelson (2002), we consider a network class with `1 regularizations on all layers except for the second last layer with an `2 regularization. With GH1 (B) = { g(x) = ReLU(vTx) : ‖v‖1 ≤ B } , a neural network class with l + 1 layers is defined as
GHl+1(B) =
{ g(x) = ReLU ( H∑ h=1 vhgh(x) ) : H∑ h=1 |vh| ≤ B, gh ∈ GHl (B) } .
Combining with the last sigmoid layer, we obtain the following discriminator class,
FHL (κ, τ,B) =
{ D(x) = sigmoid ∑ j≥1 wjsigmoid ( 2p∑ h=1 ujhgjh(x) + bj ) : ∑ j≥1 |wj | ≤ κ, 2p∑ h=1 u2jh ≤ 2, |bj | ≤ τ, gjh ∈ GHL−1(B) } .
Note that all the activation functions are ReLU(·) except that we use sigmoid(·) in the last layer in the feature map g(·). A theoretical guarantees of the class defined above is given by the following theorem.
Theorem A.1. Assume p log pn ∨ 2 ≤ c for some sufficiently small constant c > 0. Consider i.i.d. observations X1, ..., Xn ∼ (1− )N(θ, Ip) + Q and the estimator θ̂ defined by (14) with D = FHL (κ, τ,B) with H ≥ 2p, 2 ≤ L = O(1), 2 ≤ B = O(1), and τ = √ p log p. We set κ = O (√ p log p n + ) . Then, for the estimator θ̂
defined by (14) with D = FHL (κ, τ,B), we have ‖θ̂ − θ‖2 ≤ C ( p log p
n ∨ 2
) ,
with probability at least 1− e−C′(p log p+n 2) uniformly over all θ ∈ Rp such that ‖θ‖∞ ≤ √ log p and all Q.
The theorem shows that JS-GAN with a deep ReLU network can achieve the error rate p log pn ∨ 2 with respect to
the squared `2 loss. The condition ‖θ‖∞ ≤ √
log p for the ReLU network can be easily satisfied with a simple preprocessing step. We split the data into two halves, whose sizes are log n and n− log n, respectively. Then, we calculate the coordinatewise median θ̃ using the small half. It is easy to show that ‖θ̃ − θ‖∞ ≤ √ log p logn ∨
with high probability. Then, for each Xi from the second half, the conditional distribution of Xi − θ̃ given the first half is (1− )N(θ− θ̃, Ip)+ Q̃. Since √ log p logn ∨ ≤ √ log p, the condition ‖θ− θ̃‖∞ ≤ √ log p is satisfied,
and thus we can apply the estimator (14) using the shifted data Xi − θ̃ from the second half. The theoretical guarantee of Theorem A.1 will be
‖θ̂ − (θ − θ̃)‖2 ≤ C ( p log p
n ∨ 2
) ,
with high probability. Hence, we can use θ̂ + θ̃ as the final estimator to achieve the same rate in Theorem A.1.
On the other hand, our experiments show that this preprocessing step is not needed. We believe that the assumption ‖θ‖∞ ≤ √ log p is a technical artifact in the analysis of the Rademacher complexity. It can probably be dropped by a more careful analysis.

heading: B DETAILS OF EXPERIMENTS

text:
B.1 TRAINING DETAILS
The implementation for JS-GAN is given in Algorithm 1, and a simple modification of the objective function leads to that of TV-GAN. Algorithm 1 JS-GAN: argminη maxw[ 1n ∑n i=1 logDw(Xi) + E log(1−Dw(Gη(Z)))] Input: Observation set S = {X1, . . . , Xn} ∈ Rp, discriminator network Dw(x), generator network Gη(z) = z+η, learning rates γd and γg for the discriminator and the generator, batch sizem, discriminator steps in each iteration K, total epochs T , average epochs T0. Initialization: Initialize η with coordinatewise median of S. Initialize w withN(0, .05) independently on each element or Xavier (Glorot & Bengio, 2010).
1: for t = 1, . . . , T do 2: for k = 1, . . . ,K do 3: Sample mini-batch {X1, . . . , Xm} from S. Sample {Z1, . . . , Zm} from N(0, Ip) 4: gw ← ∇w[ 1mΣ m i=1 logDw(Xi) + 1 mΣ m i=1 log(1−Dw(Gη(Zi)))] 5: w ← w + γdgw 6: end for 7: Sample {Z1, . . . , Zm} from N(0, Ip) 8: gη ← ∇η[ 1mΣ m i=1 log(1−Dw(Gη(Zi)))]
9: η ← η − γggη 10: end for Return: The average estimate η over the last T0 epochs.
Several important implementation details are discussed below.
• How to tune parameters? The choice of learning rates is crucial to the convergence rate, but the minimax game is hard to evaluate. We propose a simple strategy to tune hyper-parameters including the learning rates. Suppose we have estimators θ̂1, . . . , θ̂M with corresponding discriminator networks Dŵ1 ,. . . , DŵM . Fixing η = θ̂, we further apply gradient descent to Dw with a few more epochs (but not many in order to prevent overfitting, for example 10 epochs) and select the θ̂ with the smallest value of the objective function (14) (JS-GAN) or (12) (TV-GAN). We note that training discriminator and generator alternatively usually will not suffer from overfitting since the objective function for either the discriminator or the generator is always changing. However, we must be careful about the overfitting issue when training the discriminator alone with a fixed η, and that is why we apply an early stopping strategy here. Fortunately, the experiments show if the structures of networks are same (then of course, the dimensions of the inputs are same), the choices of hyper-parameters are robust to different models and we present the critical parameters in Table 5 to reproduce the experiment results in Table 1 and Table 2.
• When to stop training? Judging convergence is a difficult task in GAN trainings, since sometimes oscillation may occur. In computer vision, people often use a task related measure and stop training once the requirement based on the measure is achieved. In our experiments below, we simply use a sufficiently large T which works well, but it is still interesting to explore an efficient early stopping rule in the future work.
• How to design the network structure? Although Theorem 3.1 and Theorem 3.2 guarantee the minimax rates of TV-GAN without hidden layer and JS-GAN with one hidden layer, one may wonder whether deeper network structures will perform better. From our preliminary experiments, TV-GAN with one hidden layer is significantly better than TV-GAN without any hidden layer. Moreover, JS-GAN

heading: 200-200-100-1 JS 0.005 0.1 2 200 25 0

text:
B.2 SETTINGS OF CONTAMINATION Q
We introduce the contamination distributions Q used in the experiments. We first consider Q = N(µ, Ip) with µ ranges in {0.2, 0.5, 1, 5}. Note that the total variation distance between N(0p, Ip) and N(µ, Ip) is of order ‖0p − µ‖ = ‖µ‖. We hope to use different levels of ‖µ‖ to test the algorithm and verify the error rate in the worst case. Second, we consider Q = N(1.5 ∗ 1p,Σ) to be a Gaussian distribution with a non-trivial covariance matrix Σ. The covariance matrix is generated according to the following steps. First generate a sparse precision matrix Γ = (γij) with each entry γij = zij ∗ τij , i ≤ j, where zij and τij are independently generated from Uniform(0.4, 0.8) and Bernoulli(0.1). We then define γij = γji for all i > j and Γ̄ = Γ + (|min eig(Γ)| + 0.05)Ip to make the precision matrix symmetric and positive definite, where min eig(Γ) is the smallest eigenvalue of Γ. The covariance matrix is Σ = Γ̄−1. Finally, we consider Q to be a Cauchy
distribution with independent component, and the jth component takes a standard Cauchy distribution with location parameter τj = 0.5.
B.3 COMPARISON DETAILS
In Section 5.3, we compare GANs with the dimension halving (Lai et al., 2016) and iterative filtering (Diakonikolas et al., 2017).
• Dimension Halving. Experiments conducted are based on the code from https://github.com/ kal2000/AgnosticMeanAndCovarianceCode. The only hyper-parameter is the threshold in the outlier removal step, and we take C = 2 as suggested in the file outRemSperical.m.
• Iterative Filtering. Experiments conducted are based on the code from https://github.com/ hoonose/robust-filter. We assume is known and take other hyper-parameters as suggested in the file filterGaussianMean.m.
B.4 SUPPLEMENTARY EXPERIMENTS FOR NETWORK STRUCTURES
The experiments are conducted with i.i.d. observations drawn from (1− )N(0p, Ip) + N(0.5 ∗ 1p, Ip) with = 0.2. Table 6 summarizes results for p = 100, n ∈ {5000, 50000} and various network structures. We observe that TV-GAN that uses neural nets with one hidden layer improves over the performance of that without any hidden layer. This indicates that the landscape of TV-GAN might be improved by a more complicated network structure. However, adding one more layer does not improve the results. For JS-GAN, we omit the results without hidden layer because of its lack of robustness (Proposition 3.1). Deeper networks sometimes improve over shallow networks, but this is not always true. We also observe that the optimal choice of the width of the hidden layer depends on the sample size.
B.5 TABLES FOR TESTING THE MINIMAX RATES
Tables 7, 8 and 9 show numerical results corresponding to Figure 2.

heading: C PROOFS OF PROPOSITION 2.1 AND PROPOSITION 3.1

text:
In the first example, consider
Q = {N(η, Ip) : η ∈ Rp} , Q̃η = {N(η̃, Ip) : ‖η̃ − η‖ ≤ r} .
In other words,Q is the class of Gaussian location family, and Q̃η is taken to be a subset in a local neighborhood of N(η, Ip). Then, with Q = N(η, Ip) and Q̃ = N(η̃, Ip), the event q̃(X)/q(X) ≥ 1 is equivalent to ‖X − η̃‖2 ≤ ‖X − η‖2. Since ‖η̃ − η‖ ≤ r, we can write η̃ = η + r̃u for some r̃ ∈ R and u ∈ Rp that satisfy 0 ≤ r̃ ≤ r and ‖u‖ = 1. Then, (8) becomes
θ̂ = argmin η∈Rp sup ‖u‖=1 0≤r̃≤r
[ 1
n n∑ i=1 I { uT (Xi − η) ≥ r̃ 2 } − P ( N(0, 1) ≥ r̃ 2 )] . (18)
Letting r → 0, we obtain (9), the exact formula of Tukey’s median.
The next example is a linear model y|X ∼ N(XT θ, 1). Consider the following classes Q = { Py,X = Py|XPX : Py|X = N(X T η, 1), η ∈ Rp } ,
Q̃η = { Py,X = Py|XPX : Py|X = N(X T η̃, 1), ‖η̃ − η‖ ≤ r } .
Here, Py,X stands for the joint distribution of y and X . The two classes Q and Q̃ share the same marginal distribution PX and the conditional distributions are specified by N(XT η, 1) and N(XT η̃, 1), respectively. Follow the same derivation of Tukey’s median, let r → 0, and we obtain the exact formula of regression depth (10). It is worth noting that the derivation of (10) does not depend on the marginal distribution PX .
The last example is on covariance/scatter matrix estimation. For this task, we set Q = {N(0,Γ) : Γ ∈ Ep}, where Ep is the class of all p×p covariance matrices. Inspired by the derivations of Tukey depth and regression depth, it is tempting to choose Q̃ in the neighborhood of N(0,Γ). However, a native choice would lead to a definition that is not even Fisher consistent. We propose a rank-one neighborhood, given by
Q̃Γ = { N(0, Γ̃) : Γ̃−1 = Γ−1 + r̃uuT ∈ Ep, |r̃| ≤ r, ‖u‖ = 1 } . (19)
Then, a direct calculation gives
I
{ dN(0, Γ̃)
dN(0,Γ) (X) ≥ 1
} = I { r̃|uTX|2 ≤ log(1 + r̃uTΓu) } . (20)
Since limr̃→0 log(1+r̃uTΓu) r̃uTΓu = 1, the limiting event of (20) is either I{|uTX|2 ≤ uTΓu} or I{|uTX|2 ≥ uTΓu}, depending on whether r̃ tends to zero from left or from right. Therefore, with the aboveQ and Q̃Γ, (8) becomes (11) under the limit r → 0. Even though the definition of (19) is given by a rank-one neighborhood of the inverse covariance matrix, the formula (11) can also be derived with Γ̃−1 = Γ−1 + r̃uuT in (19) replaced by Γ̃ = Γ + r̃uuT by applying the Sherman-Morrison formula. A similar formula to (11) in the literature is given by
Σ̂ = argmax Γ∈Ep inf ‖u‖=1
[ 1
n n∑ i=1 I{|uTXi|2 ≤ βuTΓu} ∧ 1 n n∑ i=1 I{|uTXi|2 ≥ βuTΓu}
] , (21)
which is recognized as the maximizer of what is known as the matrix depth function (Zhang, 2002; Chen et al., 2018; Paindaveine & Van Bever, 2017). The β in (21) is a scalar defined through the equation P(N(0, 1) ≤ √ β) = 3/4. It is proved in Chen et al. (2018) that Σ̂ achieves the minimax rate under Huber’s -contamination model. While the formula (11) can be derived from TV-Learning with discriminators in the form of I { dN(0,Γ̃) dN(0,Γ) (X) ≥ 1 } , a special case of (6), the formula (21) can be derived directly from TV-
GAN with discriminators in the form of I { dN(0,βΓ̃) dN(0,βΓ) (X) ≥ 1 } by following a similar rank-one neighborhood argument. This completes the derivation of Proposition 2.1.
To prove Proposition 3.1, we define F (w) = EP log sigmoid(wT g(X)) + EQ log(1 − sigmoid(wT g(X))) + log 4, so that JSg(P,Q) = maxw∈W F (w). The gradient and Hessian of F (w) are given by
∇F (w) = EP e−w T g(X)
1 + e−wT g(X) g(X)− EQ
ew T g(X)
1 + ewT g(X) g(X),
∇2F (w) = −EP ew
T g(X)
(1 + ewT g(X))2 g(X)g(X)T − EQ
e−w T g(X)
(1 + e−wT g(X))2 g(X)g(X)T .
Therefore, F (w) is concave in w, and maxw∈W F (w) is a convex optimization with a convex W . Suppose JSg(P,Q) = 0. Then maxw∈W F (w) = 0 = F (0), which implies ∇F (0) = 0, and thus we have EP g(X) = EQg(X). Now suppose EP g(X) = EQg(X), which is equivalent to ∇F (0) = 0. Therefore, w = 0 is a stationary point of a concave function, and we have JSg(P,Q) = maxw∈W F (w) = F (0) = 0.

heading: D PROOFS OF MAIN RESULTS

text:
In this section, we present proofs of all main theorems in the paper. We first establish some useful lemmas in Section D.1, and the the proofs of main theorems will be given in Section D.2.
D.1 SOME AUXILIARY LEMMAS
Lemma D.1. Given i.i.d. observations X1, ..., Xn ∼ P and the function class D defined in (13), we have for any δ > 0,
sup D∈D ∣∣∣∣∣ 1n n∑ i=1 D(Xi)− ED(X) ∣∣∣∣∣ ≤ C (√ p n + √ log(1/δ) n ) ,
with probability at least 1− δ for some universal constant C > 0. Proof. Let f(X1, ..., Xn) = supD∈D ∣∣ 1 n ∑n i=1D(Xi)− ED(X)
∣∣. It is clear that f(X1, ..., Xn) satisfies the bounded difference condition. By McDiarmid’s inequality (McDiarmid, 1989), we have
f(X1, ..., Xn) ≤ Ef(X1, ..., Xn) + √ log(1/δ)
2n ,
with probability at least 1 − δ. Using a standard symmetrization technique (Pollard, 2012), we obtain the following bound that involves Rademacher complexity,
Ef(X1, ..., Xn) ≤ 2E sup D∈D ∣∣∣∣∣ 1n n∑ i=1 iD(Xi) ∣∣∣∣∣ , (22) where 1, ..., n are independent Rademacher random variables. The Rademacher complexity can be bounded by Dudley’s integral entropy bound, which gives
E sup D∈D ∣∣∣∣∣ 1n n∑ i=1 iD(Xi) ∣∣∣∣∣ . E 1√n ∫ 2 0 √ logN (δ,D, ‖ · ‖n)dδ,
where N (δ,D, ‖ · ‖n) is the δ-covering number of D with respect to the empirical `2 distance ‖f − g‖n = √ 1 n ∑n i=1(f(Xi)− g(Xi))2. Since the VC-dimension of D is O(p), we have N (δ,D, ‖ · ‖n) . p (16e/δ)O(p) (see Theorem 2.6.7 of Van Der Vaart & Wellner (1996)). This leads to the bound 1√ n ∫ 2 0 √ logN (δ,D, ‖ · ‖n)dδ . √ p n , which gives the desired result.
Lemma D.2. Given i.i.d. observations X1, ..., Xn ∼ P, and the function class D defined in (15), we have for any δ > 0,
sup D∈D ∣∣∣∣∣ 1n n∑ i=1 logD(Xi)− E logD(X) ∣∣∣∣∣ ≤ Cκ (√ p n + √ log(1/δ) n ) ,
with probability at least 1− δ for some universal constant C > 0. Proof. Let f(X1, ..., Xn) = supD∈D ∣∣ 1 n ∑n i=1 logD(Xi)− E logD(X) ∣∣. Since sup D∈D sup x | log(2D(x))| ≤ κ,
we have sup
x1,...,xn,x′i
|f(x1, ..., xn)− f(x1, ..., xi−1, x′i, xi+1, ..., xn)| ≤ 2κ
n .
Therefore, by McDiarmid’s inequality (McDiarmid, 1989), we have f(X1, ..., Xn) ≤ Ef(X1, ..., Xn) + κ √ 2 log(1/δ)
n , (23)
with probability at least 1−δ. By the same argument of (22), it is sufficient to bound the Rademacher complexity E supD∈D ∣∣ 1 n ∑n i=1 i log(2D(Xi)) ∣∣. Since the function ψ(x) = log(2sigmoid(x)) has Lipschitz constant 1 and satisfies ψ(0) = 0, we have
E sup D∈D ∣∣∣∣∣ 1n n∑ i=1 i log(2D(Xi)) ∣∣∣∣∣ ≤ 2E sup∑ j≥1 |wj |≤κ,uj∈Rp,bj∈R ∣∣∣∣∣∣ 1n n∑ i=1 i ∑ j≥1 wjσ(u T j Xi + bj) ∣∣∣∣∣∣ , which uses Theorem 12 of Bartlett & Mendelson (2002). By Hölder’s inequality, we further have
E sup∑ j≥1 |wj |≤κ,uj∈Rp,bj∈R ∣∣∣∣∣∣ 1n n∑ i=1 i ∑ j≥1 wjσ(u T j Xi + bj) ∣∣∣∣∣∣ ≤ κEmax
j≥1 sup
uj∈Rp,bj∈R ∣∣∣∣∣ 1n n∑ i=1 iσ(u T j Xi + bj) ∣∣∣∣∣ = κE sup
u∈Rp,b∈R ∣∣∣∣∣ 1n n∑ i=1 iσ(u TXi + b) ∣∣∣∣∣ .
Note that for a monotone function σ : R→ [0, 1], the VC-dimension of the class {σ(uTx+ b) : u ∈ R, b ∈ R} is O(p). Therefore, by using the same argument of Dudley’s integral entropy bound in the proof Lemma D.1, we have
E sup u∈Rp,b∈R ∣∣∣∣∣ 1n n∑ i=1 iσ(u TXi + b) ∣∣∣∣∣ . √ p n ,
which leads to the desired result.
Lemma D.3. Given i.i.d. observations X1, .., Xn ∼ N(θ, Ip) and the function class FHL (κ, τ,B). Assume ‖θ‖∞ ≤ √ log p and set τ = √ p log p. We have for any δ > 0,
sup D∈FHL (κ,τ,B) ∣∣∣∣∣ 1n n∑ i=1 logD(Xi)− E logD(X) ∣∣∣∣∣ ≤ Cκ ( (2B)L−1 √ p log p n + √ log(1/δ) n ) ,
with probability at least 1− δ for some universal constants C > 0.
Proof. Write f(X1, ..., Xn) = supD∈FHL (κ,τ,B) ∣∣ 1 n ∑n i=1 logD(Xi)− E logD(X) ∣∣. Then, the inequality (23) holds with probability at least 1 − δ. It is sufficient to analyze the Rademacher complexity. Using the fact that the function log(2sigmoid(x)) is Lipschitz and Hölder’s inequality, we have
E sup D∈FHL (κ,τ,B) ∣∣∣∣∣ 1n n∑ i=1 i log(2D(Xi)) ∣∣∣∣∣ ≤ 2E sup
‖w‖1≤κ,‖uj∗‖2≤2,|bj |≤τ,gjh∈GHL−1(B) ∣∣∣∣∣∣ 1n n∑ i=1 i ∑ j≥1 wjsigmoid ( 2p∑ h=1 ujhgjh(Xi) + bj )∣∣∣∣∣∣ ≤ 2κE sup
‖u‖2≤2,|b|≤τ,gh∈GHL−1(B) ∣∣∣∣∣ 1n n∑ i=1 isigmoid ( 2p∑ h=1 uhgh(Xi) + b )∣∣∣∣∣ ≤ 4κE sup
‖u‖2≤2,|b|≤τ,gh∈GHL−1(B) ∣∣∣∣∣ 1n n∑ i=1 i ( 2p∑ h=1 uhgh(Xi) + b )∣∣∣∣∣ ≤ 8√pκE sup
g∈GHL−1(B) ∣∣∣∣∣ 1n n∑ i=1 ig(Xi) ∣∣∣∣∣+ 4κτE ∣∣∣∣∣ 1n n∑ i=1 i ∣∣∣∣∣ . Now we use the notation Zi = Xi − θ ∼ N(0, Ip) for i = 1, ..., n. We bound E supg∈GHL−1(B) ∣∣ 1 n ∑n i=1 ig(Zi + θ)
∣∣ by induction. Since E ( sup
g∈GH1 (B)
1
n n∑ i=1 ig(Zi + θ)
)
≤ E ( sup ‖v‖1≤B 1 n n∑ i=1 iv T (Zi + θ) )
≤ B ( E ∣∣∣∣∣ 1n n∑ i=1 iZi ∣∣∣∣∣ ∞ + ‖θ‖∞E ∣∣∣∣∣ 1n n∑ i=1 i ∣∣∣∣∣ ) ≤ CB √
log p+ ‖θ‖∞√ n ,
and
E ( sup
g∈GHl+1(B)
1
n n∑ i=1 ig(Zi + θ)
)
≤ E ( sup
‖v‖1≤B,gh∈GHl (B)
1
n n∑ i=1 i H∑ h=1 vhgh(Zi + θ)
)
≤ BE ( sup
g∈GHl (B) ∣∣∣∣∣ 1n n∑ i=1 ig(Zi + θ) ∣∣∣∣∣ )
≤ 2BE ( sup
g∈GHl (B)
1
n n∑ i=1 ig(Zi + θ)
) ,
we have
E ( sup
g∈GHL−1(B)
1
n n∑ i=1 ig(Zi + θ)
) ≤ C(2B)L−1 √ log p+ ‖θ‖∞√
n .
Combining the above inequalities, we get
E ( sup
D∈FHL (κ,τ,B)
1
n n∑ i=1 i logD(Zi + θ)
) ≤ Cκ ( √ p(2B)L−1 √ log p+ ‖θ‖∞√
n + τ√ n
) .
This leads to the desired result under the conditions on τ and ‖θ‖∞.
D.2 PROOFS OF MAIN THEOREMS
Proof of Theorem 3.1. We first introduce some notations. Define F (P, η) = maxw,b Fw,b(P, η), where
Fw,b(P, η) = EP sigmoid(w TX + b)− EN(η,Ip)sigmoid(w TX + b).
With this definition, we have θ̂ = argminη F (Pn, η), where we use Pn for the empirical distribution 1 n ∑n i=1 δXi . We shorthand N(η, Ip) by Pη , and then
F (Pθ, θ̂) ≤ F ((1− )Pθ + Q, θ̂) + (24)
≤ F (Pn, θ̂) + + C
(√ p
n +
√ log(1/δ)
n
) (25)
≤ F (Pn, θ) + + C
(√ p
n +
√ log(1/δ)
n
) (26)
≤ F ((1− )Pθ + Q, θ) + + 2C
(√ p
n +
√ log(1/δ)
n
) (27)
≤ F (Pθ, θ) + 2 + 2C
(√ p
n +
√ log(1/δ)
n
) (28)
= 2 + 2C
(√ p
n +
√ log(1/δ)
n
) . (29)
With probability at least 1− δ, the above inequalities hold. We will explain each inequality. Since
F ((1− )Pθ + Q, η) = max w,b [(1− )Fw,b(Pθ, η) + Fw,b(Q, η)] ,
we have sup η |F ((1− )Pθ + Q, η)− F (Pθ, η)| ≤ ,
which implies (24) and (28). The inequalities (25) and (27) are implied by Lemma D.1 and the fact that
sup η |F (Pn, η)− F ((1− )Pθ + Q, η)| ≤ sup w,b ∣∣∣∣∣ 1n n∑ i=1 sigmoid(wTXi + b)− Esigmoid(wTX + b) ∣∣∣∣∣ . The inequality (26) is a direct consequence of the definition of θ̂. Finally, it is easy to see that F (Pθ, θ) = 0, which gives (29). In summary, we have derived that with probability at least 1− δ,
Fw,b(Pθ, θ̂) ≤ 2 + 2C
(√ p
n +
√ log(1/δ)
n
) ,
for all w ∈ Rp and b ∈ R. For any u ∈ Rp such that ‖u‖ = 1, we take w = u and b = −uT θ, and we have
f(0)− f(uT (θ − θ̂)) ≤ 2 + 2C
(√ p
n +
√ log(1/δ)
n
) ,
where f(t) = ∫
1 1+ez+tφ(z)dz, with φ(·) being the probability density function of N(0, 1). It is not hard to
see that as long as |f(t)− f(0)| ≤ c for some sufficiently small constant c > 0, then |f(t)− f(0)| ≥ c′|t| for some constant c′ > 0. This implies
‖θ̂ − θ‖ = sup ‖u‖=1 |uT (θ̂ − θ)|
≤ 1 c′ sup ‖u‖=1 ∣∣∣f(0)− f(uT (θ − θ̂))∣∣∣ . + √ p
n +
√ log(1/δ)
n ,
with probability at least 1− δ. The proof is complete.
Proof of Theorem 3.2. We continue to use Pη to denote N(η, Ip). Define
F (P, η) = max ‖w‖1≤κ,u,b Fw,u,b(P, η),
where Fw,u,b(P, η) = EP logD(X) + EN(η,Ip) log (1−D(X)) + log 4,
with D(x) = sigmoid (∑
j≥1 wjσ(u T j x+ bj)
) . Then,
F (Pθ, θ̂) ≤ F ((1− )Pθ + Q, θ̂) + 2κ (30)
≤ F (Pn, θ̂) + 2κ + Cκ
(√ p
n +
√ log(1/δ)
n
) (31)
≤ F (Pn, θ) + 2κ + Cκ
(√ p
n +
√ log(1/δ)
n
) (32)
≤ F ((1− )Pθ + Q, θ) + 2κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) (33)
≤ F (Pθ, θ) + 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) (34)
= 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) .
The inequalities (30)-(34) follow similar arguments for (24)-(28). To be specific, (31) and (33) are implied by Lemma D.2, and (32) is a direct consequence of the definition of θ̂. To see (30) and (34), note that for any w such that ‖w‖1 ≤ κ, we have
| log(2D(X))| ≤ ∣∣∣∣∣∣ ∑ j≥1 wjσ(u T j X + bj) ∣∣∣∣∣∣ ≤ κ. A similar argument gives the same bound for | log(2(1−D(X)))|. This leads to
sup η |F ((1− )Pθ + Q, η)− F (Pθ, η)| ≤ 2κ , (35)
which further implies (30) and (34). To summarize, we have derived that with probability at least 1− δ,
Fw,u,b(Pθ, θ̂) ≤ 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) ,
for all ‖w‖1 ≤ κ, ‖uj‖ ≤ 1 and bj . Take w1 = κ, wj = 0 for all j > 1, u1 = u for some unit vector u and b1 = −uT θ, and we get
fuT (θ̂−θ)(κ) ≤ 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) , (36)
where fδ(t) = E log
2
1 + e−tσ(Z) + E log
2
1 + etσ(Z+δ) , (37)
with Z ∼ N(0, 1). Direct calculations give
f ′δ(t) = E e−tσ(Z) 1 + e−tσ(Z) σ(Z)− E e tσ(Z+δ) 1 + etσ(Z+δ) σ(Z + δ), f ′′δ (t) = −Eσ(Z)2 e−tσ(Z)
(1 + e−tσ(Z))2 − Eσ(Z + δ)2 e
tσ(Z+δ)
(1 + etσ(Z+δ))2 . (38)
Therefore, fδ(0) = 0, f ′δ(0) = 1 2 (Eσ(Z)− Eσ(Z + δ)), and f ′′ δ (t) ≥ − 12 . By the inequality
fδ(κ) ≥ fδ(0) + κf ′δ(0)− 1
4 κ2,
we have κf ′δ(0) ≤ fδ(κ) + κ2/4. In view of (36), we have
κ
2
(∫ σ(z)φ(z)dz − ∫ σ(z + uT (θ̂ − θ))φ(z)dz ) ≤ 4κ + 2Cκ (√ p
n +
√ log(1/δ)
n
) + κ2
4 .
It is easy to see that for the choices of σ(·), ∫ σ(z)φ(z)dz − ∫ σ(z + t)φ(z)dz is locally linear with respect to t. This implies that
κ‖θ̂ − θ‖ = κ sup ‖u‖=1 uT (θ̂ − θ) . κ
( + √ p
n +
√ log(1/δ)
n
) + κ2.
Therefore, with a κ . √
p n + , the proof is complete.
Proof of Theorem 4.1. We use Pθ,Σ,h to denote the elliptical distribution EC(θ,Σ, h). Define
F (P, (η,Γ, g)) = max ‖w‖1≤κ,u,b Fw,u,b(P, (η,Γ, g)),
where Fw,u,b(P, (η,Γ, g)) = EP logD(X) + EEC(η,Γ,g) log (1−D(X)) + log 4,
with D(x) = sigmoid (∑
j≥1 wjσ(u T j x+ bj)
) . Let P be the data generating process that satisfies
TV(P, Pθ,Σ,h) ≤ , and then there exist probability distributions Q1 and Q2, such that
P + Q1 = Pθ,Σ,h + Q2.
The explicit construction of Q1, Q2 is given in the proof of Theorem 5.1 of Chen et al. (2018). This implies that
|F (P, (η,Γ, g))− F (Pθ,Σ,h, (η,Γ, g))| ≤ sup
‖w‖1≤κ,u,b |Fw,u,b(P, (η,Γ, g))− Fw,u,b(Pθ,Σ,h, (η,Γ, g))|
= sup ‖w‖1≤κ,u,b |EQ2 log(2D(X))− EQ1 log(2D(X))| ≤ 2κ . (39)
Then, the same argument in Theorem 3.2 (with (35) replaced by (39)) leads to the fact that with probability at least 1− δ,
Fw,u,b(Pθ,Σ,h, (θ̂, Σ̂, ĥ)) ≤ 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) ,
for all ‖w‖1 ≤ κ, ‖uj‖ ≤ 1 and bj . Take w1 = κ, wj = 0 for all j > 1, u1 = u/ √ uT Σ̂u for some unit vector
u and b1 = −uT θ/ √ uT Σ̂u, and we get
fuT (θ̂−θ)√ uT Σ̂u
(κ) ≤ 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) ,
where
fδ(t) =
∫ log ( 2
1 + e−tσ(∆s)
) h(s)ds+ ∫ log ( 2
1 + etσ(δ+s)
) ĥ(s)ds,
where δ = u T (θ̂−θ)√ uT Σ̂u and ∆ = √ uTΣu√ uT Σ̂u . A similar argument to the proof of Theorem 3.2 gives
κ
2
(∫ σ(∆s)h(s)ds− ∫ σ(δ + s)ĥ(s)ds ) ≤ 4κ + 2Cκ (√ p
n +
√ log(1/δ)
n
) + κ2
4 .
Since ∫ σ(∆s)h(s)ds = 1
2 =
∫ σ(s)ĥ(s)ds,
the above bound is equivalent to
κ 2 (H(0)−H(δ)) ≤ 4κ + 2Cκ
(√ p
n +
√ log(1/δ)
n
) + κ2
4 ,
where H(δ) = ∫ σ(δ+ s)ĥ(s)ds. The above bound also holds for κ2 (H(δ)−H(0)) by a symmetric argument, and therefore the same bound holds for κ2 |H(δ)−H(0)|. Since H ′(0) = ∫ σ(s)(1− σ(s))ĥ(s)ds = 1, H(δ) is locally linear at δ = 0, which leads to a desired bound for δ = u T (θ̂−θ)√ uT Σ̂u . Finally, since uT Σ̂u ≤ M , we get the bound for uT (θ̂ − θ). The proof is complete by taking supreme over all unit vector u.
Proof of Theorem A.1. We continue to use Pη to denote N(η, Ip). Define
F (P, η) = sup D∈FHL (κ,τ,B) FD(P, η),
with FD(P, η) = EP logD(X) + EN(η,Ip) log(1−D(X)) + log 4.
Follow the same argument in the proof of Theorem 3.2, use Lemma D.3, and we have
FD(Pθ, θ̂) ≤ Cκ ( + (2B)L−1 √ p log p
n +
√ log(1/δ)
n
) ,
uniformly over D ∈ FHL (κ, τ,B) with probability at least 1 − δ. Choose w1 = κ and wj = 0 for all wj > 1. For any unit vector ũ ∈ Rp, take u1h = −u1(h+p) = ũh for h = 1, ..., p and b1 = −ũT θ. For h = 1, ..., p, set g1h(x) = max(xh, 0). For h = p + 1, ..., 2p, set g1h(x) = max(−xh−p, 0). It is obvious that such u and b satisfy ∑ h u 2 1h ≤ 2 and |b1| ≤ ‖θ‖ ≤ √ p‖θ‖∞ ≤ √ p log p. We need to show both the functions max(x, 0) and max(−x, 0) are elements of GHL−1(B). This can be proved by induction. It is obvious that max(xh, 0),max(−xh, 0) ∈ GH1 (B) for any h = 1, ..., p. Suppose we have max(xh, 0),max(−xh, 0) ∈ GHl (B) for any h = 1, ..., p. Then,
max (max(xh, 0)−max(−xh, 0), 0) = max(xh, 0), max (max(−xh, 0)−max(xh, 0), 0) = max(−xh, 0).
Therefore, max(xh, 0),max(−xh, 0) ∈ GHl+1(B) as long as B ≥ 2. Hence, the above construction satisfies D(x) = sigmoid(κsigmoid(ũT (x− θ))) ∈ FHL (κ, τ,B), and we have
fuT (θ̂−θ)(κ) ≤ Cκ
( + (2B)L−1 √ p log p
n +
√ log(1/δ)
n
) , (40)
where the definition of fδ(t) is given by (37) with Z ∼ N(0, 1) and σ(·) is taken as sigmoid(·). Apply the a similar in the proof of Theorem 3.2, we obtain the desired result.


summaries: In this paper, the authors proposed a new format of representation called scene programs, to describe the visual scenes. To extract the scene programs from scenes, the authors exploited the off-the-shelf object detection and segmentations model, mask r-cnn to extract all objects and the corresponding attributes from the images, and then detect groups for those objects, which are then used to generate the programs which matches the input scenes. The experiments are performed on a synthetics datasets which consists of multiple shapes with different attributes. The experiments shows that the proposed model can infer more accurate programs from the scenes, and those generated programs can be used to recover the input scenes more accurately. Besides, the authors also showed that the generated scene programs can be used for image editing and making visual analogy.

summaries: This paper presents a system that infers programs describing 3D scenes composed of simple primitives. The system consists of three stages each of which is trained separately. First, the perceptual module extracts object masks and their attributes. The objects are then are split into several groups. Finally, each group is mapped to a corresponding DSL program using a sequence-to-sequence network similar to the ones typically employed in neural machine translation.

abstractText: Human scene perception goes beyond recognizing a collection of objects and their pairwise relations. We understand higher-level, abstract regularities within the scene such as symmetry and repetition. Current vision recognition modules and scene representations fall short in this dimension. In this paper, we present scene programs, representing a scene via a symbolic program for its objects, attributes, and their relations. We also propose a model that infers such scene programs by exploiting a hierarchical, object-based scene representation. Experiments demonstrate that our model works well on synthetic data and transfers to real images with such compositional structure. The use of scene programs has enabled a number of applications, such as complex visual analogy-making and scene extrapolation.

authors:
Jimmy Ba
Volodymyr Mnih
Koray Kavukcuoglu

title: Multiple object recognition with visual attention

venue: In ICLR,

authors:
Dzmitry Bahdanau
Kyunghyun Cho
Yoshua Bengio

title: Neural machine translation by jointly learning to align and translate

title: Pix2code: Generating code from a graphical user interface screenshot

venue: In ACM SIGCHI Symposium on Engineering Interactive Computing Systems,

authors:
Rudy Bunel
Matthew Hausknecht
Jacob Devlin
Rishabh Singh
Pushmeet Kohli

title: Leveraging grammar and reinforcement learning for neural program synthesis

authors:
Yuntian Deng
Anssi Kanervisto
Jeffrey Ling
Alexander M Rush

title: Image-to-markup generation with coarse-to-fine attention

venue: In ICML,

authors:
Jacob Devlin
Jonathan Uesato
Surya Bhupatiraju
Rishabh Singh
Abdel-rahman Mohamed
Pushmeet Kohli

title: Robustfill: Neural program learning under noisy i/o

authors:
Kevin Ellis
Daniel Ritchie
Armando Solar-Lezama
Josh Tenenbaum

title: Learning to infer graphics programs from hand-drawn images

venue: In NeurIPS,

authors:
SM Eslami
Nicolas Heess
Theophane Weber
Yuval Tassa
Koray Kavukcuoglu
Geoffrey E Hinton

title: Attend, infer, repeat: Fast scene understanding with generative models

venue: NeurIPS,

authors:
Yaroslav Ganin
Tejas Kulkarni
Igor Babuschkin
S.M. Ali Eslami
Oriol Vinyals

title: Synthesizing programs for images using reinforced adversarial learning

authors:
Kaiming He
Xiangyu Zhang
Shaoqing Ren
Jian Sun

title: Deep residual learning for image recognition

venue: In CVPR,

authors:
Geoffrey E Hinton
Ruslan R Salakhutdinov

title: Reducing the dimensionality of data with neural networks

authors:
Jonathan Huang
Kevin Murphy

title: Efficient inference in occlusion-aware generative models of images

venue: In ICLR Workshop,

authors:
Phillip Isola
Jun-Yan Zhu
Tinghui Zhou
Alexei A Efros

title: Image-to-image translation with conditional adversarial networks

authors:
Justin Johnson
Ranjay Krishna
Michael Stark
Li-Jia Li
David Shamma
Michael Bernstein
Li Fei-Fei

title: Image retrieval using scene graphs

authors:
Justin Johnson
Bharath Hariharan
Laurens van der Maaten
Li Fei-Fei
C Lawrence Zitnick
Ross Girshick

title: Clevr: A diagnostic dataset for compositional language and elementary visual reasoning

authors:
Thomas N Kipf
Ethan Fetaya
Kuan-Chieh Wang
Max Welling
Richard S Zemel

title: Neural relational inference for interacting systems

authors:
Jun Li
Kai Xu
Siddhartha Chaudhuri
Ersin Yumer
Hao Zhang
Leonidas Guibas

title: GRASS: Generative Recursive Autoencoders for Shape Structures

venue: In SIGGRAPH,

authors:
Manyi Li
Akshay Gadi Patil
Kai Xu
Siddhartha Chaudhuri
Owais Khan
Ariel Shamir
Changhe Tu
Baoquan Chen
Daniel Cohen-Or
Hao Zhang

title: Grains: Generative recursive autoencoders for indoor scenes

venue: ACM TOG,

authors:
Thang Luong
Hieu Pham
Christopher D. Manning

title: Effective approaches to attention-based neural machine translation

venue: In EMNLP,

authors:
Chengjie Niu
Jun Li
Kai Xu

title: Im2Struct: Recovering 3D Shape Structure from a Single RGB Image

authors:
Emilio Parisotto
Abdel-rahman Mohamed
Rishabh Singh
Lihong Li
Dengyong Zhou
Pushmeet Kohli

title: Neuro-symbolic program synthesis

authors:
Scott E Reed
Yi Zhang
Yuting Zhang
Honglak Lee

title: Deep visual analogy-making

authors:
Irvin Rock
Stephen Palmer

title: The legacy of gestalt psychology

venue: Csgnet:

heading: 1 INTRODUCTION

text:
When examining the image in Figure 1a, we instantly recognize the shape, color, and material of the objects it depicts. We can also effortlessly imagine how we may extrapolate the set of objects in the scene while preserving object patterns (Figure 1b). Our ability to imagine unseen objects arises from holistic scene perception: we not only recognize individual objects from an image, but naturally perceive how they should be organized into higher-level structure (Rock & Palmer, 1990).
Recent AI systems for scene understanding have made impressive progress on detecting, segmenting, and recognizing individual objects (He et al., 2017). In contrast, the problem of understanding high-level, abstract relations among objects is less studied. While a few recent papers have attempted to produce a holistic scene representation for scenes with a variable number of objects (Ba et al., 2015; Huang & Murphy, 2015; Eslami et al., 2016; Wu et al., 2017), the relationships among these objects are not captured in these models.
The idea of jointly discovering objects and their relations has been explored only very recently, where the learned relations are often in the form of interaction graphs (van Steenkiste et al., 2018; Kipf et al., 2018) or semantic scene graphs (Johnson et al., 2015), both restricted to pairwise, local relations. However, our ability to imagine extrapolated images as in Figure 1 relies on our knowledge of long-range, hierarchical relationships among objects, such as how objects are grouped and what patterns characterize those groups.
In this paper, we aim to tackle the problem of understanding higher-level, abstract regularities such as repetition and symmetry. We propose to represent scenes as scene programs. We define a domainspecific language for scenes, capturing both objects with their geometric and semantic attributes, as well as program commands such as loops to enforce higher-level structural relationships. Given an image of a complex scene, we propose to infer its scene program via a hierarchical bottom-up approach. First, we parse the image into individual objects and infer their attributes, resulting in the object representation. Then, we organize these objects into different groups, i.e. the group representation, where objects in each group fall into the same program block. Finally, we describe each group with a program, and combine these programs to get the program representation for the entire scene.
∗This work was done when Yunchao Liu was a visiting student at MIT CSAIL.
Our model applies deep neural networks for each stage of this process and is able to generate programs describing the input image with high accuracy. When testing on scenes that are more complex than those used for training, our hierarchical inference process achieves better generalization performance than baseline methods that attempt to infer a program directly from the image. Our model is also able to handle ambiguity, generating multiple possible programs when there is more than one way to describe the scene. Furthermore, our method generalizes to real-world images without any additional supervised training programs; only the low-level object detection module must be re-trained. Finally, we demonstrate how our model facilitates high-level image editing, as users can change parameters in the inferred program to achieve the editing effects they want more efficiently. We show examples of such image edits, including extrapolations such as the one in (Figure 1b), on both synthetic and photographic images.
Our contributions are therefore three-fold:
1. We propose scene programs: a new representation for scenes, drawing insights from classic findings in cognitive science and computer graphics.
2. We present a method for inferring scene programs from images using a hierarchical approach (from objects to groups to programs).
3. We demonstrate that our model can achieve high accuracy on describing both synthetic and constrained real scenes with programs. Combined with modern image-to-image translation methods, our model generates realistic images of extrapolated scenes, capturing both highlevel scene structure and low-level object appearance.

heading: 2 RELATED WORK

text:
Describing Images with Programs Ellis et al. (2018) performs a similar task as ours where handdrawn images of 2D geometry primitives are converted to high-level programs. This work uses a constraint-based SAT solver to perform program search and is much slower than neural network models. IM2LATEX (Deng et al., 2017) de-renders images into low-level LATEX markup using a neural network, while our work discovers high-level programs from an image of objects. SPIRAL (Ganin et al., 2018) uses reinforcement learning to infer a sequence of low-level drawing commands that can reproduce an image, meanwhile learning a distribution from which images can be sampled. Beltramelli (2018) learns to convert GUI images to markup-like code. Unlike these papers, our model performs program induction in 3D and infers high-level structural patterns both in object layout and color.
Describing the Structure of 3D Shapes and Scenes Beyond 2D images, prior work in vision and graphics has attempted to infer high-level structure from 3D objects and 3D scenes. The most relevant to our approach are those that extract a so-called symmetry hierarchy, in which 3D geometry is hierarchically grouped by either attachment or symmetric relationships (Wang et al., 2011). This
representation has been used to train generative models of 3D shapes (Li et al., 2017) and indoor 3D scenes (Li et al., 2019), as well as to infer a hierarchical bounding box structure from a single image of a 3D shape (Niu et al., 2018). Our program representation bears some resemblance to the symmetry hierarchy, but it generalizes to repetitive patterns beyond symmetries and also models patterns in object visual attributes (e.g., color). CSGNet (Sharma et al., 2018) learns to parse shapes with a set of primitive and arithmetic commands; Tulsiani et al. (2017) parses shapes into an assembly of geometric primitives. In this paper, we focus on learning the high-level scene regularities described by loop structures.
Neural Program Synthesis In general, a program synthesis model outputs an explicit program by learning from examples. Recent works on neural program synthesis include R3NN (Parisotto et al., 2017) and RobustFill (Devlin et al., 2017), which perform end-to-end program synthesis from input/output examples. Bunel et al. (2018) goes beyond the pure supervised learning setting and improves performance on diversity and syntax by leveraging grammar and reinforcement learning. These models synthesize programs based on input/output pairs, which is different from our setting, where a program is generated to describe an input image. In the vision domain, Sun et al. (2018) learns decision strategies represented as programs from demonstration videos, while we focus on describing the complex correlations among objects in static scenes.

heading: 3 METHOD

text: Our model combines vision and sequence models via structured representations. An object parser predicts the segmentation mask and attributes for each object in the image. A group recognizer predicts the group that each object belongs to. Finally, a program synthesizer generates a program block for each object group. Figure 2 shows an example of synthesizing programs from an input image, where a sphere is selected at random (highlighted) and the group that this object belongs to is predicted, which consists of six spheres. Then the program for this group (highlighted) is synthesized.

heading: 3.1 A DOMAIN-SPECIFIC LANGUAGE (DSL) FOR SCENES

text:
In order to constrain the program space to make it tractable for our models, we introduce human prior on scene regularities that can be described as programs. More specifically, we introduce a Domain Specific Language (DSL) which explicitly defines the space of our scene programs. We present the grammar of our DSL in Table 1, which contains 3 primitive commands (cube, sphere, cylinder) and 2 loop structures (for, rotate). The positions for each object are defined as affine transformations of loop indices, while the colors are more complicated functions of the loop indices, displaying alternating (modular) and repeating (division) patterns.
Furthermore, since the DSL allows unbounded program depth, we define program blocks to further reduce complexity. Each type of program block is an production instance of the Statement token, and
objects that belong to the same block form a group. For example, in this work the program blocks include single objects, layered for loops of depth ≤ 3, and single-layer rotations of ≤ 4 objects.

heading: 3.2 OBJECT PARSING

text: Following the spirit of The Trace Hypothesis (Ellis et al., 2018), we use object attributes as an intermediate representation between image space and structured program space. Parsing individual objects from the input image consists of two steps: mask prediction and attribute prediction. For each object, its instance segmentation mask is predicted by a Mask R-CNN (He et al., 2017). Next, the mask is concatenated with the original image, and sent to a ResNet-34 (He et al., 2015) to predict object attributes. In our work, object attributes include shape, size, material, color and 3D coordinates. Each attribute is encoded as a one-hot vector, except for coordinates. The overall representation of an object is a vector of length 18. The networks are trained with ground truth masks and attributes, respectively. For the attribute network, we minimize the mean-squared error between output and ground truth attributes.

heading: 3.3 GROUP DETECTION

text:
When we identify a distinct visual pattern, we first know which objects in the image form the pattern before we can tell what the pattern is. Motivated by this idea, we develop a group recognizer that tells us which objects form a group that can be described by a single program block. The group recognizer works after mask prediction is performed, and answers the following specific question: given an input object, which objects are in the same group with this object?
The input to the model consists of three parts: the original image, the mask of the input object, and the mask of all objects. These three parts are concatenated and sent to a ResNet-152 followed by fully connected layers. The output contains two parts: a binary vector g where g[i] = 1 denotes object i in the same group with the input object, and the category c of the group, representing the type of program block that this group belongs to. The network is trained to minimize the binary cross entropy loss for group recognition, and the cross entropy loss for category classification.

heading: 3.4 NEURAL PROGRAM SYNTHESIS

text:
With the object attributes and groups obtained from the vision models, the final step in our model is to generate program sequences describing the input image. Since we have already detected object groups, what remains is to generate a program block for each group. For this goal we train a sequence to sequence (seq2seq) LSTM with an encoder-decoder structure and attention mechanism (Luong et al., 2015; Bahdanau et al., 2015). The input sequence is a set of object attributes that form a group, which are sorted by their 3D coordinates. The output program consists of two parts: program tokens are predicted as a sequence as in neural machine translation, and program parameters are predicted by a MLP from the hidden state at each time step. At each step, we predict a token t as well as a
Algorithm 1: Combining group prediction with program synthesis Result: a program sequence P Input: a set of object attributes O; while O is not empty do
randomly choose oi ∈ O; predict the group that contains oi, indexed by G; also predict the group category c; get attributes of objects that belong to the group, A = {oj |j ∈ G}; remove A from O; send A, c to program synthesizer, get program p; add p to P ;
end
parameter matrix P , which contains predicted parameters for all possible tokens. Then we use P [t] as the output parameter for this step.
Since the program synthesizer only works for a single group, a method for combining the group prediction with program synthesis is needed. Consider the simplest case where we randomly choose an object and describe the group it belongs to. This procedure is described in Algorithm 1. In practice, by default we sample 10 times and stop when a correct program is generated. Here correct means that we can recover the scene attributes successfully by executing the program.

heading: 4 EXPERIMENTS

text: We perform several experiments on synthetic scene images, including quantitative comparison with baseline methods and further extensions and applications. We further demonstrate our model’s ability to generalize to real images with a small amount of hand-labeled supervision which is only at the object level. We also apply our method to other tasks, specifically image extrapolation and visual analogy-making, on both synthetic and real images.

heading: 4.1 DATASET

text:
We create a synthetic dataset of images rendered from complex scenes with rich program structures. Figure 3 displays some examples drawn from the dataset. These images are generated by first sampling scenes and then rendering using the same renderer as in CLEVR (Johnson et al., 2017). Each scene consists of a few groups, where objects in the same group can be described by a program block. The groups are sampled from predefined program primitives with multi-layered translational and rotational symmetries. Further, we also incorporate rich color patterns into the primitives.
Our synthetic dataset includes annotations on object attributes and programs. We train and test the models on two synthetic datasets, REGULAR and RANDOM, each containing 20,000 training and 500 test images, where each image has at most 2 groups of multiple objects, in addition to many groups of a single object. In the REGULAR dataset (Figure 3a), the objects are placed on a grid and have discrete coordinates. We increase the scene complexity by adding randomness in the RANDOM dataset (Figure 3c), where objects are placed uniformly at random with continuous coordinates and different sizes.

heading: 4.2 VISUAL PROGRAM SYNTHESIS

text:
Setup. We present evaluation results on our synthetic dataset introduced above. We compare with an ablated version of our full model, where we use a simple search-based heuristic grouping method (HG) which does not require any training or any knowledge of the program patterns. The details of this method are presented in Appendix A. We also compare with two baselines. One removes group recognition and instead synthesizes programs from all object attributes (derender-LSTM). Another directly synthesizes programs from the input image in an end-to-end manner (CNN-LSTM). The model uses a CNN as encoder and a LSTM with attention as decoder. We use the same network
architecture as in attention-based neural image captioning (Xu et al., 2015), except that the decoder predicts a token as well as a parameter matrix at each time step.
We evaluate the models on both the REGULAR and the RANDOM test sets. For evaluation on generalization, we also create an additional test set of 100 images, where each image contains three groups of multiple objects. These images are more complex and harder to describe than those in training.
Results. Figure 3 includes qualitative results generated by our model on all three test sets. Our model generates accurate results in the REGULAR setting and is able to recognize the two groups from neighbouring objects (Figure 3a). Although trained on images with two groups, our model can perform well when tested on images with three groups (Figure 3b). When objects are placed at random, our model can accurately recognize which objects form a regular pattern and describe them with programs (Figure 3c). For quantitative evaluation, we compute program token accuracy and parameter loss, defined as the percentage of correctly predicted tokens and the mean-squared error of parameter prediction, respectively. To evaluate the global performance of the generated program, we also compute reconstruction accuracy of the programs, defined as the percentage of programs that correctly reconstruct the original image. The reconstruction accuracy is evaluated on all three test sets.
We present the test results in Table 2, where our model outperforms baseline methods in each of the metrics, and achieves good performance on generalization. Note that the deep grouping model outperforms the simple heuristic grouping method, as it learns from the data distribution specified by
our program space. Also note that our model performs better on RANDOM, while the baselines do not perform as well. This is because our group detection model discovers groups among randomly placed objects better than among regularly placed objects, as it is easier to rule out outsiders when they look random.
Tackling ambiguous input. While our model can generate program representations for images with high accuracy, it can also generate multiple possible programs when the input is ambiguous. Figure 4 shows an example where the red group can be described by either a two-layer for loop or a rotation of 4 objects. Our hierarchical method allows explicit specification of group category. When executing Algorithm 1, instead of selecting the most confident group category, we search top 3 proposals, and execute the synthesized program block to decide if each proposal corresponds to a possible correct program. Figure 4 demonstrates programs generated by our model, while the baseline methods tend to collapse to one possible answer and is unable to generate others.
Program synthesis from partial observations. Our model can also handle scenes where there are invisible (or hardly visible) objects. Figure 5 demonstrates how our model operates on these scenes. Given an input image, we generate object instance masks and remove those with area below a certain threshold, so that the remaining objects can be correctly recognized. These objects form the partial observation of our model, from which the program synthesizer generates a program block which correctly describes the scene, including (partially) occluded objects. The flexibility of the neural program synthesizer allows us to recognize the same program pattern given different partial observations. Consider the two examples at the bottom of Figure 5. They have different set of observations (8 and 6 objects on bottom left and right, respectively) due to the different distances, and our model is able to correctly recognize both of them.

heading: 4.3 IMAGE EDITING

text:
Image editing via program representation. With the expressive power of the program representation, our model can be applied to tasks that require a high-level structural knowledge of the scene. For example, when an image lies within the space defined by our DSL, it can be efficiently edited using the program representation generated by our model. Figure 6 shows some examples of image editing, where the input image (Figure 6a) is represented by a program. Users can then edit the program to achieve the preferred editing effects. The edited program is sent to a graphics engine to render the new image. The structural form of our program representation allows various types of high-level editing, including spacial extrapolation (Figure 6b, c), changing color patterns (Figure 6d) and shapes (Figure 6e). Each of the four examples requires only one edit in the program, while using the traditional object representation, users have to change objects one at a time, averaging 6.25 edits per image.
Real image extrapolation. An advantage of our method which uses object attributes as a connection between vision and program synthesis is to generalize to real images. Since our neural program synthesizer is independent from visual recognition, only the vision systems need to be retrained for our entire model to work on real images.
Figure 7a shows images of LEGO blocks shot from a camera in real-world settings. We create a dataset of 120 real images, where we use 90 for training, 10 for validation, and 20 for testing. To adapt our model to generate programs for these images, we first pretrain on a synthetic dataset of 4,000 images rendered by a graphics engine. Then we fine-tune the model on 90 real images with labeled masks and attributes. The vision system is then linked with the pretrained program synthesizer which does not require any fine-tuning. Even with a small amount of real data for fine-tuning, our model generalizes well and correctly predicts the programs for each test image.
Furthermore, the image editing techniques introduced above can also be applied to such real images. Here we present an experiment on real image extrapolation. Given an input image, we generate the program describing the image and also extract object patches with Mask R-CNN. The program is extended by increasing the iteration number, which is a simple way of “imagining” what could be the next given a sequence of observations.
Our original method uses a graphics engine to render new images from edited programs (Figure 6), which is not applicable for real images. For this purpose, we use pix2pix (Isola et al., 2017) as an approximate neural renderer. After program inference, we execute the edited program and retrieve newly added object masks. These masks can be computed using camera parameters and 3D coordinates, while here we use retrieval for simplicity. All of the patches are pasted on a white background, and then sent to pix2pix to generate realistic background and lighting. Figure 7c displays the editing results. The edited images preserve object appearances in the original images, and also fix the errors made by mask prediction (small white gaps in Figure 7b) and contain realistic-looking shadows.

heading: 4.4 VISUAL ANALOGY MAKING

text:
Besides representing images for efficient editing, scene programs can also be used as encoded images. For example, the distance in program space can also be applied to model similarity between images, which is already introduced by (Ellis et al., 2018). Motivated by this idea, we consider visual analogy making (Reed et al., 2015), where an input image is converted to a new image given other reference images. We introduce a setting where the reference is an image pair and ask the intuitive question, if B follows A, then what should follow C?
Here we use a simple solution based on representation distance. More specifically, for an encoder R and an input image c with reference pair (a, b), we set R(d) = R(c)+R(b)−R(a) and decode R(d) to get the output. In our case, the encoder is our program synthesis model, while we use pix2pix as a neural decoder. In order to perform arithmetic operations, the program is represented as a matrix, where each line starts with a token followed by parameters (see Appendix A.3 for details). We compare our model with an autoencoder (Hinton & Salakhutdinov, 2006). The autoencoder we adopt takes an input image of size 256× 256, encodes the input into a 256-dimensional vector and then decodes the encoded vector back to original image size.
Figure 8 shows qualitative results of the visual analogy making task. Using our program representation, our model generates perceptually plausible results (Figure 8e). While the autoencoder
can sometimes correctly change the number of objects, it fails to preserve the layout arrangements (Figure 8f). We also compute the average L2 distance between model output and ground truth images made by humans. Table 3 shows that our model generates images that are closer to the ground truth than the baseline.

heading: 5 CONCLUSION

text:
We propose scene programs as a structured representation of complex scenes with high-level regularities. We also present a novel method that infers scene programs from 2D images in a hierarchical bottom-up manner. Our model achieves high accuracy on a synthetic dataset and also generalizes to real images. The representation power of programs allows our model to be applied to other tasks in computer vision, such as image editing and analogy making, on both synthetic and photographic images.
Acknowledgements. We thank Jiayuan Mao for insightful discussions and anonymous reviewers for their helpful feedback. This work was supported in part by NSF #1231216, NSF #1447476, NSF #1753684, ONR MURI N00014-16-1-2007, Facebook, and the Yao Class Exchange Program at IIIS, Tsinghua University.