Robust Learning of Fixed-Structure Bayesian Networks

Introduction

Probabilistic graphical models [KF09] provide an appealing and unifying formalism to succinctly represent structured high-dimensional distributions. The general problem of inference in graphical models is of fundamental importance and arises in many applications across several scientific disciplines, see [WJ08] and references therein. In this work, we study the problem of learning graphical models from data [Nea03, DSA11]. There are several variants of this general learning problem depending on: (i) the precise family of graphical models considered (e.g., directed, undirected), (ii) whether the data is fully or partially observable, and (iii) whether the structure of the underlying graph is known a priori or not (parameter estimation versus structure learning). This learning problem has been studied extensively along these axes during the past five decades, (see, e.g., [CL68, Das97, AKN06, WRL06, AHHK12, SW12, LW12, BMS13, BGS14, Bre15]) resulting in a beautiful theory and a collection of algorithms in various settings.

The main vulnerability of all these algorithmic techniques is that they crucially rely on the assumption that the samples are precisely generated by a graphical model in the given family. This simplifying assumption is inherent for known guarantees in the following sense: if there exists even a very small fraction of arbitrary outliers in the dataset, the performance of known algorithms can be totally compromised. It is important to explore the natural setting when the aforementioned assumption holds only in an approximate sense. Specifically, we study the following broad question:

Can we efficiently learn graphical models when a constant fraction of the samples are corrupted, or equivalently, when the model is slightly misspecified?

In this paper, we focus on the model of corruptions considered in [DKK+16] (Definition 1) which generalizes many other existing models, including Huber’s contamination model [Hub64]. Intuitively, given a set of good samples (from the true model), an adversary is allowed to inspect the samples before corrupting them, both by adding corrupted points and deleting good samples. In contrast, in Huber’s model, the adversary is oblivious to the samples and is only allowed to add bad points.

We would like to design robust learning algorithms for Question 1 whose sample complexity, $N$, is close to the information-theoretic minimum, and whose computational complexity is polynomial in $N$. We emphasize that the crucial requirement is that the error guarantee of the algorithm is independent of the dimensionality $d$ of the problem.

In this work, we study Question 1 in the context of Bayesian networks [JN07]. We focus on the fully observable case when the underlying network is given. In the non-robust setting, this learning problem is straightforward: the “empirical estimator” (which coincides with the maximum likelihood estimator) is known to be sample and computationally efficient [Das97]. In sharp contrast, even this most basic regime is surprisingly challenging in the robust setting. For example, the very special case of robustly learning a Bernoulli product distribution (corresponding to an empty network with no edges) was analyzed only recently in [DKK+16].

To formally state our results, we first give a detailed description of the corruption model we study.

Given $0<\epsilon<1/2$ and a distribution family $\mathcal{P}$, the algorithm specifies some number of samples $N$, and $N$ samples $X_{1},X_{2},\ldots,X_{N}$ are drawn from some (unknown) ground-truth $P\in\mathcal{P}$. The adversary is allowed to inspect $P$ and the samples, and replaces $\epsilon N$ of them with arbitrary points. The set of $N$ points is then given to the algorithm. We say that a set of samples is $\epsilon$-corrupted if it is generated by this process.

Note that $P$ is determined by $G$ and $p$. We will frequently index $p$ as a vector. We use the notation $p_{k}$ and the associated events $\Pi_{k}$, where each $k\in[m]$ stands for an $(i,a)\in\Gamma$ lexicographically ordered.

Our Results.

We give the first efficient robust learning algorithm for Bayesian networks with a known graph $G$. Our algorithm has information-theoretically near-optimal sample complexity, runs in time polynomial in the size of the input (the samples), and provides an error guarantee that scales near-linearly with the fraction of adversarially corrupted samples, under the following restrictions: First, we assume that each parental configuration is reasonably likely. Intuitively, this assumption seems necessary because we need to observe each configuration many times in order to learn the associated conditional probability to good accuracy. Second, we assume that each of the conditional probabilities is balanced, i.e., bounded away from and $1$. This assumption is needed for technical reasons. In particular, we need this to show that a good approximation to the conditional probability table implies that the corresponding Bayesian network is close in total variation distance.

Our algorithm is given in Section 3. We first note that the sample complexity of our algorithm is near-optimal for learning Bayesian networks with known structure. The following sample complexity lower bound holds even without corrupted samples:

Let $\mathcal{BN}_{d,f}$ denote the family of Bernoulli Bayesian networks on $d$ variables such that every node has at most $f$ parents. The worst-case sample complexity of learning $\mathcal{BN}_{d,f}$, within total variation distance $\epsilon$ and with probability $9/10$, is $\Omega(2^{f}\cdot d/\epsilon^{2})$ for all $f\leq d/2$ when the graph structure is known.

Consider Bayes nets whose average in-degree is close to the maximum in-degree, that is, when $m=\Theta(2^{f}d)$, the sample complexity lower bound in Fact 4 becomes $\Omega(m/\epsilon^{2})$, so our sample complexity is optimal up to polylogarithmic factors.

We remark that Theorem 3 is most useful when $c$ is a constant and the Bayesian network has bounded fan-in $f$. In this case, the condition on $\alpha$ follows from the $c$-balanced assumption: When both $c$ and $f$ are constants, $\alpha=\min_{(i,a)\in S}\Pr_{P}[\Pi_{i,a}]\geq c^{f}$ is also a constant, so the condition $c^{f}\geq\Omega(\epsilon\sqrt{\log(1/\epsilon)})$ automatically hold when $\epsilon$ is smaller than some constant. On the other hand, the problem of learning Bayesian networks is less interesting when the fan-in is too large. For example, if some node has $f=\omega(\log(d))$ parents, then the size of the conditional probability table is at least $2^{f}$, which is super-polynomial in the dimension $d$.

Experiments.

We performed an experimental evaluation of our algorithm on both synthetic and real data. Our evaluation allowed us to verify the accuracy and the sample complexity rates of our theoretical results. In all cases, the experiments validate the usefulness of our algorithm, which significantly outperforms previous approaches, almost exactly matching the best rate without noise.

Related Work.

Question 1 fits in the framework of robust statistics [HR09, HRRS86]. Classical estimators from this field can be classified into two categories: either (i) they are computationally efficient but incur an error that scales polynomially with the dimension $d$, or (ii) they are provably robust (in the aforementioned sense) but are hard to compute. In particular, essentially all known estimators in robust statistics (e.g., the Tukey median [Tuk75]) have been shown [JP78, Ber06, HM13] to be intractable in the high-dimensional setting. We note that the robustness requirement does not typically pose information-theoretic impediments for the learning problem. In most cases of interest (see, e.g., [CGR15, CGR16, DKK+16]), the sample complexity of robust learning is comparable to its (easier) non-robust variant. The challenge is to design computationally efficient algorithms.

Efficient robust estimators are known for various low-dimensional structured distributions (see, e.g., [DDS14, CDSS13, CDSS14a, CDSS14b, ADLS16, ADLS17, DLS18]). However, the robust learning problem becomes surprisingly challenging in high dimensions. Recently, there has been algorithmic progress on this front: [DKK+16, LRV16] give polynomial-time algorithms with improved error guarantees for certain “simple” high-dimensional structured distributions. The results of [DKK+16] apply to simple distributions, including Bernoulli product distributions, Gaussians, and mixtures thereof (under some natural restrictions). Since the works of [DKK+16, LRV16], computationally efficient robust estimation in high dimensions has received considerable attention (see, e.g., [DKS17, DKK+17, BDLS17, DKK+18a, DKS18b, DKS18a, HL18, KSS18, PSBR18, DKK+18b, KKM18, DKS18c, LSLC18]).

2 Overview of Algorithmic Techniques

Our algorithmic approach builds on the framework of [DKK+16] with new technical and conceptual ideas. At a high level, our algorithm works as follows: We draw an $\epsilon$-corrupted set of samples from a Bayesian network $P$ with known structure, and then iteratively remove samples until we can return the empirical conditional probability table.

First, we associate a vector $F(X)$ to each sample $X$ so that learning the mean of $F(X)$ to good accuracy is sufficient to recover the distribution. In the case of binary products, $F(X)$ is simply $X$, while in our case we need to take into account additional information about conditional means.

From this point, our algorithm will try to do one of two things: Either we show that the sample mean of $F(X)$ is close to the conditional mean of the true distribution (in which case we can already learn the ground-truth Bayes net $P$), or we are able to produce a filter, i.e., we can remove some of our samples, and it is guaranteed that we throw away more bad samples than good ones. If we produce a filter, we then iterate on those samples that pass the filter. To produce a filter, we compute a matrix $M$ which is roughly the empirical covariance matrix of $F(X)$. We show that if the corruptions are sufficient to notably disrupt the sample mean of $F(X)$, there must be many erroneous samples that are all far from the mean in roughly the same direction, and we can detect this direction by looking at the largest eigenvector of $M$. If we project all samples onto this direction, concentration bounds of $F(X)$ will imply that almost all samples far from the mean are erroneous, and thus filtering them out will provide a cleaner set of samples.

Section 2 contains some technical results specific to Bayesian networks that we need. Section 3 gives the details of our algorithm and an overview of its analysis. In Section 4, we present the experimental evaluations. In Section 5, we conclude and propose directions for future work.

Technical Preliminaries

The structure of this section is as follows: First, we bound the total variation distance between two Bayes nets in terms of their conditional probability tables. Second, we define a function $F(x,q)$, which takes a sample $x$ and returns an $m$-dimensional vector that contains information about the conditional means. Finally, we derive a concentration bound from Azuma’s inequality. Proofs from this section have been deferred to Appendix A.

Lemma 5 says that to learn a balanced fixed-structure Bayesian network, it is sufficient to learn all the relevant conditional means. However, each sample $x\sim P$ gives us information about $p_{i,a}$ only if $x\in\Pi_{i,a}$. To resolve this, we map each sample $x$ to an $m$-dimensional vector $F(x,q)$, and “fill in” the entries that correspond to conditional means for which the condition failed to happen. We will set these coordinates to their empirical conditional means $q$:

When $q=p$ (the true conditional means), the expectation of the $(i,a)$-th coordinate of $F(X,p)$, for $X\sim P$, is the same conditioned on either $\Pi_{i,a}$ or $\neg\Pi_{i,a}$. Using the conditional independence properties of Bayesian networks, we will show that the covariance of $F(x,p)$ is diagonal.

Finally, we will need a suitable concentration inequality that works under conditional independence properties. We can use Azuma’s inequality to show that the projections of $F(X,q)$ on any direction $v$ is concentrated around the projection of the sample mean $q$.

Robust Learning Algorithm

We first look into the major ingredients required for our filtering algorithm, and compare our proof with that for product distributions in [DKK+16] on a more technical level.

In Section 2, we mapped each sample $X$ to $F(X,q)$ which contains information about the conditional means $q$, and we showed that it is sufficient to learn the mean of $F(X,q)$ to learn the ground-truth Bayes net.

Let $M$ denote the empirical covariance matrix of $(F(X,q)-q)$. We decompose $M$ into three parts: One coming from the ground-truth distribution, one coming from the subtractive error (because the adversary can remove $\epsilon N$ good samples), and one coming from the additive error (because the adversary can add $\epsilon N$ bad samples). We will make use of the following observations:

The noise-free distribution has a diagonal covariance matrix.

The term coming from the subtractive error has no large eigenvalues.

These two observations imply that any large eigenvalues of $M$ are due to the additive error. Finally, we will reuse our concentration bounds to show that if the additive errors are frequently far from the mean in a known direction, then they can be reliably distinguished from good samples.

For the case of binary product distributions in [DKK+16], (1) is trivial because the coordinates are independent; but for Bayesian networks we need to expand the dimension of the samples and fill in the missing entries properly. Condition (2) is due to concentration bounds, and for product distributions it follows from standard Chernoff bounds, while for Bayes nets, we must instead rely on martingale arguments and Azuma’s inequality. The main difference between the proof of correctness of our algorithm and those given in [DKK+16] lies in analyzing the mean, covariance and tail bounds of $F(X,q)$, and showing that its mean and covariance are well-behaved when $q$ is close to $p$ (see Lemma 24 in Appendix B).

First, we need to show that a large enough set of samples with no noise satisfy properties we expect from a representative set of samples. We need that the mean, covariance, and tail bounds of $F(X,p)$ behave like we would expect them to. This happens with high probability. The details are given in Lemma 24 in Appendix B.1. We call a set of samples that satisfies these properties $\epsilon$-good for $P$.

Our algorithm takes as input an $\epsilon$-corrupted multiset $S^{\prime}$ of $N={\widetilde{\Omega}}(m\log(1/\tau)/\epsilon^{2})$ samples. We write $S^{\prime}=(S\setminus L)\cup E$, where $S$ is the set of samples before corruption, $L$ contains the good samples that have been removed or (in later iterations) incorrectly rejected by filters, and $E$ represents the remaining corrupted samples. We assume that $S$ is $\epsilon$-good. In the beginning, we have $|E|+|L|\leq 2\epsilon|S|$. As we add filters in each iteration, $E$ gets smaller and $L$ gets larger. However, we will prove that our filter rejects more samples from $E$ than $S$, so $|E|+|L|$ must get smaller.

We will prove Theorem 3 by iteratively running the following efficient filtering procedure:

Let $0<\epsilon<1/2$. Let $P$ be a $c$-balanced Bayesian network on $\{0,1\}^{d}$ with known structure $G$. Assume each parental configuration of $P$ occurs with probability at least $\alpha\geq\Omega(\epsilon\sqrt{\log(1/\epsilon)}/c)$. Let $S^{\prime}=S\cup E\setminus L$ be a set of samples such that $S$ is $\epsilon$-good for $P$ and $|E|+|L|\leq 2\epsilon|S^{\prime}|$. There is an algorithm that, given $G$, $\epsilon$, and $S^{\prime}$, runs in time ${\widetilde{O}}(d|S^{\prime}|)$, and either

Returns an $S^{\prime\prime}=S\cup E^{\prime}\setminus L^{\prime}$ such that $|S^{\prime\prime}|\leq(1-\frac{\epsilon}{d\ln d})|S^{\prime}|$ and $|E^{\prime}|+|L^{\prime}|<|E|+|L|$.

If this algorithm produces a subset $S^{\prime\prime}$, then we iterate using $S^{\prime\prime}$ in place of $S^{\prime}$. We will present the algorithm establishing Proposition 9 in the following section. We first use it to prove Theorem 3.

Next we analyze the running time. Observe that we can filter out at most $2\epsilon N$ samples, because we reject more bad samples than good ones. By Proposition 9, every time we produce a filter, we remove at least ${\widetilde{\Omega}}(d/\epsilon)|S^{\prime}|={\widetilde{\Omega}}(Nd/\epsilon)$ samples. Therefore, there are at most ${\widetilde{O}}(d)$ iterations, and each iteration takes time ${\widetilde{O}}(d|S^{\prime}|)={\widetilde{O}}(Nd)$ by Proposition 9, so the overall running time is ${\widetilde{O}}(Nd^{2})$. ∎

2 Algorithm Filter-Known-Topology

In this section, we present Algorithm 1 that establishes Proposition 9. We use $X\in_{u}S$ to denote that the point $X$ is drawn uniformly from the set of samples $S$.

At a high level, Algorithm 1 computes a matrix $M$, and shows that: either $\|M\|_{2}$ is small, and we can output the empirical conditional probabilities, or $\|M\|_{2}$ is large, and we can use the top eigenvector of $M$ to remove bad samples.

Our first step is to analyze the spectrum of $M$, and in particular show that $M$ is close in spectral norm to $w_{E}M_{E}$. To do this, we begin by showing that the spectral norm of $M_{S,0}$ is relatively small. Since $S$ is good, we have bounds on the second moments $F(X,p)$. We just need to deal with the error from replacing $p$ with $q$ (see Appendix B.2 for the proof):

$\|M_{S,0}\|_{2}\leq O(\epsilon+\sqrt{\sum_{k}\Pr_{S}[\Pi_{k}](p_{k}-q_{k})^{2}}+\sum_{k}\Pr_{S}[\Pi_{k}](p_{k}-q_{k})^{2})$.

Next, we wish to bound the contribution to $M$ coming from the subtractive error. We show that this is small due to concentration bounds on $P$ and hence on $S$. The idea is that for any unit vector $v$, we have tail bounds for the random variable $v\cdot(F(X,q)-q)$ and, since $L$ is a subset of $S$, $L$ can at worst consist of a small fraction of the tail of this distribution.

$w_{L}\|M_{L}\|_{2}\leq O(\epsilon\log(1/\epsilon)+\epsilon\|p-q\|_{2}^{2})$.

Finally, combining the above results, since $M_{S}$ and $M_{L}$ have small contribution to the spectral norm of $M$ when $\|p-q\|_{2}$ is small, most of it must come from $M_{E}$.

$\|M-w_{E}M_{E}\|_{2}\leq O\left(\epsilon\log(1/\epsilon)+\sqrt{\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2}}+\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2}\right)$.

Lemma 12 follows using the identity $|S^{\prime}|M=|S|M_{S,0}+|E|M_{E,0}-|L|M_{L,0}$ and bounding the errors due to the diagonals of $M_{E}$ and $M_{L}$.

The Case of Small Spectral Norm.

$\sqrt{\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}]^{2}(p_{k}-q_{k})^{2}}\leq 2\sqrt{\epsilon\|M\|_{2}}+O(\epsilon\sqrt{\log(1/\epsilon)+1/\alpha})$.

If $\|M\|_{2}\leq O(\epsilon\log(1/\epsilon)/\alpha)$, then

The Case of Large Spectral Norm.

Now we consider the case when $\|M\|_{2}\geq C\epsilon\ln(1/\epsilon)/\alpha$. We begin by showing that $p$ and $q$ are not too far apart from each other. The bound given by Lemma 15 is now dominated by the $\|M\|_{2}$ term. Lower bounding the $\Pr_{S^{\prime}}[\Pi_{k}]$ by $\alpha$ gives the following claim.

$\|p-q\|_{2}\leq\delta:=3\sqrt{\epsilon\|M\|_{2}}/\alpha$.

Recall that $v^{*}$ is the largest eigenvector of $M$. We project all the points $F(X,q)$ onto the direction of $v^{*}$. Next we show that most of the variance of $(v^{*}\cdot(F(X,q)-q))$ comes from $E$.

$v^{\ast T}(w_{E}M_{E})v^{\ast}\geq\frac{1}{2}v^{\ast T}Mv^{\ast}$.

Claim 18 follows from the observation that $\|M-w_{E}M_{E}\|_{2}$ is much smaller than $\|M\|_{2}$. This is obtained by substituting the bound on $\|p-q\|_{2}$ (in terms of $\|M\|_{2}$) from Claim 17 into the bound on $\|M-w_{E}M_{E}\|_{2}$ given by Lemma 12.

Claim 18 implies that the tails of $w_{E}E$ are reasonably thick. In particular, we show that there must be a threshold $T>0$ satisfying the desired property in Step 9 of our algorithm.

If Lemma 19 were not true, by integrating this tail bound, we can show that $v^{\ast T}M_{E}v^{\ast}$ would be small. Therefore, Step 11 of Algorithm 11 is guaranteed to find some valid threshold $T>0$.

Finally, we show that the set of samples $S^{\prime\prime}$ we return after the filter is better than $S^{\prime}$ in terms of $|L|+|E|$. This completes the proof of the second case of Proposition 9.

If we write $S^{\prime\prime}=S\cup E^{\prime}\setminus L^{\prime}$, then $|E^{\prime}|+|L^{\prime}|<|E|+|L|$ and $|S^{\prime\prime}|\leq(1-\frac{\epsilon}{d\ln d})|S^{\prime}|$.

Claim 9 follows from the fact that $S$ is $\epsilon$-good, so we only remove at most $(3\exp(T^{2}/2)+\epsilon/T^{2}\log d)|S|$ samples from $S$. Since we remove more than twice as many samples from $S^{\prime}$, most of the samples we throw away are from $E$. Moreover, we remove at least $(1-\frac{\epsilon}{d\ln d})|S^{\prime}|$ samples because we can show that the threshold $T$ is at most $\sqrt{d}$.

Running Time of Our Algorithm 1

Experiments

We test our algorithms using data generated from both synthetic and real-world networks (e.g., the ALARM network [BSCC89]) with synthetic noise. All experiments were run on a laptop with 2.6 GHz CPU and 8 GB of RAM. We found that our algorithm achieves the smallest error consistently in all trials, and that the error of our algorithm almost matches the error of the empirical conditional probabilities of the uncorrupted samples. Moreover, our algorithm can easily scale to thousands of dimensions with millions of samples. The bottleneck of our algorithm is fitting millions of samples of thousands dimension all in the memory.

We draw the parameters of $P$ independently from $[0,1/4]\cup[3/4,1]$ uniformly at random, i.e., in a setting where the “balancedness” assumption does not hold. Our experiments show that our filtering algorithm works very well in this setting, even when the assumptions under which we can prove theoretical guarantees are not satisfied. This complements our theoretical results and illustrates that our algorithm is not limited by these assumptions and can apply to more general settings in practice.

In Figure 1, we compare the performance of (1) our filtering algorithm, (2) the empirical conditional probability table with noise, and (3) a RANSAC-based algorithm (see the end of Section 4 for a detailed description). We use the error of the empirical conditional mean without noise (i.e., MLE estimator with only good samples) as the gold standard, since this is the best one could hope for even if all the corrupted samples are identified. We tried various graph structures for the Bayes net $P$ and noise distributions, and similar patterns arise for all of them. In the top figure, the dependency graph of $P$ is a randomly generated tree, and the noise distribution is a binary product distribution; In the bottom figure, the dependency graph of $P$ is a random graph, and the noise distribution is the tree Bayes net used as the ground truth in the first experiment. The reader is referred to Appendix C.1 for a full description of how we generate the dependency graphs and noise distributions.

2 Semi-Synthetic Experiments

In the semi-synthetic experiments, we apply our algorithm to robustly learn real-world Bayesian networks. The ALARM network [BSCC89] is a classic Bayes net that implements a medical diagnostic system for patient monitoring.

Our experimental setup is as follows: The underlying graph of ALARM has $37$ nodes and $509$ parameters. Since the variables in ALARM can have up to $4$ different values, we first transform it into an equivalent binary-valued Bayes net(see Appendix C.3 for more details). After the transformation, the network has $d=61$ nodes and $m=820$ parameters. We are interested in whether our filtering algorithm can learn a Bayes net that is “close” to ALARM when samples are corrupted; and how many corrupted samples can our algorithm tolerate. For $\epsilon=[0.05,0.1,\ldots,0.4]$, we draw $N=10^{6}$ samples, where a $(1-\epsilon)$-fraction of the samples come from ALARM, and the other $\epsilon$-fraction comes from a noise distribution.

In Figure 2, we compare the performance of (1) our filtering algorithm, (2) the empirical conditional means with noise, and (3) a RANSAC-based algorithm. We use the error of the empirical conditional means without noise as the gold standard. We tried various noise distributions and observed similar patterns. In Figure 2, the noise distribution is a Bayes net with random dependency graphs and conditional probabilities drawn from $[0,\frac{1}{4}]\cup[\frac{3}{4},1]$ (same as the ground-truth Bayes net in Figure 1).

The experiments show that our filtering algorithm outperforms MLE and RANSAC, and that the error of our algorithm degrades gracefully as $\epsilon$ increases. It is worth noting that even the ALARM network does not satisfy our balancedness assumption on the parameters, our algorithm still performs well on it and recovers the conditional probability table of ALARM in the presence of corrupted samples.

RANSAC uses subsampling in the hope of getting an estimator that is not affected too much by the noise. The hope is that a small subsample might not contain any erroneous points. The approach proceeds by computing many such estimators and appropriately selecting the best one.

In our experiments, we let RANSAC select $10\%$ of the samples uniformly at random, and repeat this process $100$ times. After a subset of samples are selected, we compute the empirical conditional means and estimate the total variation distance between the corresponding Bayes net and the ground truth. Since we know the ground truth, we can make it easier for RANSAC by selecting the best hypothesis that it ever produced during its execution.

The main conceptual message of our experimental evaluation of RANSAC is that it does not perform well in high dimensions for the following reason: To guarantee that there are very few noisy points for such a subsample, we must take an exponential (in the dimension) number of subsets. We are not the first to observe that RANSAC does not work for robustly learning high-dimensional distributions. Previously, [DKK+17] showed that RANSAC does not work in practice for the problem of robustly learning a spherical Gaussian.

Conclusions and Future Directions

In this paper, we initiated the study of the efficient robust learning for graphical models. We described a computationally efficient algorithm for robustly learning Bayesian networks with a known topology, under some mild assumptions on the conditional probability table. We evaluate our algorithm experimentally, and we view our experiments as a proof of concept demonstration that our techniques can be practical for learning fixed-structure Bayesian networks. A challenging open problem is to generalize our results to the case when the underlying directed graph is unknown.

This work is part of a broader agenda of systematically investigating the robust learnability of high-dimensional structured probability distributions. There is a wealth of natural probabilistic models that merit investigation in the robust setting, including undirected graphical models (e.g., Ising models), and graphical models with hidden variables (i.e., incorporating latent structure).

We are grateful to Daniel Hsu for suggesting the model of Bayes nets, and for pointing us to [Das97]. Yu Cheng is supported in part by NSF CCF-1527084, CCF-1535972, CCF-1637397, CCF-1704656, IIS-1447554, and NSF CAREER Award CCF-1750140. Ilias Diakonikolas is supported by NSF CAREER Award CCF-1652862 and a Sloan Research Fellowship. Daniel Kane is supported by NSF CAREER Award CCF-1553288 and a Sloan Research Fellowship.

References

Appendix A Omitted Proofs from Section 2

In this section, we give proofs for the technical lemmas in Section 2. Lemma 5 bounds the total variation distance between two balanced Bayesian networks in terms of their conditional probability tables. Lemma 5 is a simple corollary of Lemma 21.

Let $P$ and $Q$ be Bayesian networks with the same dependency graph $G$. In terms of the conditional probability tables $p$ and $q$ of $P$ and $Q$, we have:

Let $A$ and $B$ be two distributions on $\{0,1\}^{d}$. We have:

Fix $i\in[d]$. The events $\Pi_{i,a}$ form a disjoint partition of $\{0,1\}^{d}$. Dividing the sum above into this partition, we obtain

Let $P_{\leq i}$ and $Q_{\leq i}$ be the distribution over the first $i$ coordinates of $P$ and $Q$ respectively. Let $P_{i}$ and $Q_{i}$ be the distribution of the $i$-th coordinate of $P$ and $Q$ respectively.

The first and the fifth steps use Equation 2, the second step uses that the $i$-th coordinate is independent of the first $(i-1)$ coordinates conditioned on $\Pi_{i,a}$, and the third and fourth steps use Equation 1.

Now observe that the $P_{i}\mid\Pi_{i,a}$ and $Q_{i}\mid\Pi_{i,a}$ are Bernoulli distributions with means $p_{i,a}$ and $q_{i,a}$. For $p,q\in$, we have:

Lemma 5 gives a simpler expression for total variation distance between two $c$-balanced binary Bayesian networks whose minimum probability of any $\Pi_{k}$ is at least $\epsilon$.

We associate a vector $F(X)$ to each sample $X$, so that $F(X)$ contains information about the conditional means, and learning the mean of $F(X)$ to good accuracy is sufficient to recover the distribution. Recall that $q$ is the vector of empirical conditional means, and we define $F(x,q):\{0,1\}^{d}\to^{m}$ as follows (Definition 6): If $x\in\Pi_{i,a}$, then $F(x,q)_{i,a}=x_{i}$, otherwise $F(x,q)_{i,a}=q_{i,a}$.

We will prove some properties of $F$. First, we note that $F$ is invertible in the following sense.

Fix $q\in^{m}$ and $j\in[d]$. Given $(x_{1},\ldots,x_{j})$, we can compute $F(x,q)_{i,a}$ for all $(i,a)$ with $i\leq j$. We can recover $(x_{1},\dots,x_{j})$ from these $F(x,q)_{i,a}$ as well.

By the definition of $F(x,q)$, to compute $F(x,q)_{i,a}$ we need to know $x_{i}$ and whether $x\in\Pi_{i,a}$. Note that whether or not $x\in\Pi_{i,a}$ depends only on $(x_{1},\dots,x_{i-1})$, so $F(x,p)_{i,a}$ is a function of $(x_{1},\dots,x_{i-1},x_{i})$.

We will show by induction that $(x_{1},\dots,x_{j})$ can be recovered from all $F(x,p)_{i,a}$ with $i\leq j$. Since $x_{1}$ has no parents, we have $x_{1}=F(x,p)_{1,a^{\prime}}$ for the empty bitstring $a^{\prime}$. For $i>1$, we have that $x_{i}=F(x,p)_{i,a}$ for the unique $a$ with $x\in\Pi_{i,a}$, and we can decide which $a$ based on $(x_{1},\dots,x_{i-1})$. ∎

Next, we show that when $q=p$ (the true conditional probabilities), although the coordinates of $F(X,p)$ are not independent, the mean of a coordinate of $F(X,p)$ remains unchanged even if we condition on the values of previous coordinates.

Let $k=(i,a)$. Since we order the $(i,a)$’s lexicographically, $(F(X,p)_{1},\dots,F(X,p)_{k-1})$ includes $F(X)_{j,a^{\prime}}$ for all $(j,a^{\prime})$ with $j<i$. By Claim 22, these determine the value of the parents of $X_{i}$, i.e., whether or not $\Pi_{i,a}$ occurs.

We build on Claim 23 to show that, although the coordinates of $F(X,p)$ are not independent, the first and second moments are the same as that of a product distribution of the marginal of each coordinate.

Finally, we will need a suitable concentration inequality that works under conditional independence. Lemma 8 shows that the projections of $F(X,q)$ on any direction $v$ is concentrated around its mean.

Consider an $x\in\{0,1\}^{d}$. If $x\in\Pi_{i,a}$, then we have $F(x,p)_{i,a}=F(x,q)_{i,a}=x_{i}$ and so

If $x\notin\Pi_{i,a}$, then $F(x,p)_{i,a}=p_{i,a}$ and $F(x,q)_{i,a}=q_{i,a}$, hence

An application of the Cauchy-Schwarz inequality gives that, if $|v\cdot(F(X,q)-q)|\geq T+\|p-q\|_{2}$ then $|v\cdot(F(x,p)-p)|\geq T$. Therefore, the probability of the former holding for $X$ must be at most the probability that the latter holds for $X$. ∎

Appendix B Omitted Proofs from Sections 3

This section analyzes Algorithm 1 and gives the proof of Proposition 9.

The basic idea of the analysis is as follows: If the empirical conditional probability table $q$ is close to the true conditional probability table $p$ of $P$, then outputting $q$ is correct. We know that we have enough samples that the empirical conditional probability table with no noise ${\widetilde{p}}$ is a good approximation to $p$. Therefore, we will be in good shape so long as the corruption of our samples does not introduce a large error in the conditional probability table.

Thinking more concretely about this error, we may split it into two parts: $L$, the subtractive error, and $E$ the additive error. Using concentration results for $P$, it can be shown that the subtractive errors cannot cause significant problems for the conditional probability table. It remains to consider additive errors. The bad samples in $E$ can introduce notable errors in the conditional probability table, since any given sample can be $\sqrt{d}$ far from the mean. If many of the corrupted samples line up in the same direction, this can lead to a notable discrepancy.

However, if many of these errors line up in some direction (which is necessary in order to have a large impact on the mean), the effects will be reflected in the first two moments. More concretely, if for some unit vector $v$, the expectation of $v\cdot F(E,q)$ is very far from the expectation of $v\cdot F(P,q)$, this will force the variance of $v\cdot F(S^{\prime},q)$ to be large. This implies two things: First, it tells us that if $v\cdot F(S^{\prime},q)$ is small for all $v$ (a condition equivalent to $\|M\|_{2}$ being small), we know that $q$ is a good approximation to the true conditional probability table. Second, if $\|M\|_{2}$ is large, we can find a unit vector $v$ where $v\cdot F(S^{\prime},q)$ has large variance. A reasonable fraction of this variance must be coming from samples in $E$ that have $v\cdot F(X,q)$ very far from the mean. On the other hand, using concentration bounds for $v\cdot F(P,q)$, we know that very few valid samples are this far from the mean. This discrepancy will allow us to create a filter which rejects more samples from $E$ than from $S$.

In Section B.1, we will provide a set of deterministic conditions that we expect from the good samples and show that they happen with high probability. In Section B.2, we will prove some structural lemmas about the spectrum of $M$. In Section B.3, we will show that if $\|M\|_{2}$ is small, then we can output the empirical conditional probabilities. In Section B.4, we will show that if $\|M\|_{2}$ is large, then we can use the top eigenvector of $M$ to remove bad samples.

Given a large enough set $S$ of good samples drawn from the ground-truth Bayesian network $P$, Lemma 24 states that, for $X\in_{u}S$, the mean, covariance, and tail bounds of $F(X,p)$ behave like we would expect them to. We call a set of samples that satisfies these properties $\epsilon$-good for $P$.

Let $S$ be a set of $\Omega((m\log(m/\epsilon)+\log(1/\tau))\cdot\log^{2}d\cdot\epsilon^{-2})$ samples from $P$. Let $p$ and ${\widetilde{p}}$ denote the conditional probability tables of $P$ and of the empirical distribution given by $S$ respectively. Then, with probability at least $1-\tau$, we have the following:

$|\Pr_{S}[\Pi_{k}]-\Pr_{P}[\Pi_{k}]|\leq\epsilon$,

$\sum_{k}\Pr_{S}[\Pi_{k}]({\widetilde{p}}_{k}-p_{k})^{2}\leq\epsilon^{2}$,

For all unit vectors $v$ and $T>0$, we have

For (i), by the Chernoff and union bounds, with probability at least $1-\tau/10$, we have that our empirical estimates for $\Pr_{P}[\Pi_{k}]$ are correct to within $\epsilon$ as long as we have at least $O(\log(m/\tau)/\epsilon^{2})$ samples.

Note that for a fixed $k=(i,a)$, we have ${\widetilde{p}}_{k}$ is the empirical expectation of $\Pr_{S}[\Pi_{k}]N$ independent samples from a Bernoulli with probability $p_{k}$. By Chernoff bounds, when $N\geq\Omega(m\log(m/\tau)/\epsilon^{2})$, we have $|{\widetilde{p}}_{k}-p_{k}|\leq\epsilon/\sqrt{m\Pr_{S}[\Pi_{k}]}$ with probability at least $1-\tau/10m$, . By a union bound, this holds for all $k$ except with probability at most $1/10\tau$. Then we have $\sum_{k}\Pr_{S}[\Pi_{k}]({\widetilde{p}}_{k}-p_{k})^{2}\leq\epsilon^{2}$.

For (iii) and (iv), we first prove this happens for a fixed $v$ and $T$ with sufficiently high probability and then take a union bound over a cover of $v$ and $T$.

Let $S$ be a set of $N$ samples from $P$. Let $X\in_{u}S$ and $Y\sim P$. For any unit vector $v$ and $T\geq 0$, we have that

$\Pr[|v\cdot(F(X,p)-p)|\geq T]\leq 5\exp(-T^{2}/2)/2+\epsilon/(2T^{2})$,

with probability at least $1-\exp(-\Omega(N\epsilon^{2}))$.

By Lemma 8, we have $\Pr[|v\cdot(F(Y,p)-p)|\geq T]\leq 2\exp(-T^{2}/2)$. Hence,

is the sum of $N$ i.i.d. Bernoulli random variables, each with mean at most $2\exp(-T^{2}/2)$. We use the following two versions of the Chernoff bound:

Let $Z_{1},\ldots,Z_{N}$ be i.i.d. Bernoullis with mean $\mu$. Then

$\Pr[\sum_{i}Z_{i}/N\geq(1+\delta)\mu]\leq\exp(-\delta\ln(1+\delta)N\mu/2)$ for $\delta>0$.

$\Pr[\sum_{i}Z_{i}/N\geq\nu]\leq\exp(-D(\nu||\mu)N)$ for $\nu\geq\mu$, where $D(\nu||\mu)=\nu\ln(\nu/\mu)+(1-\nu)\ln((1-\nu)/(1-\mu))$ is the KL-divergence between Bernoullis with probabilities $\nu$ and $\mu$ .

Here $Z_{i}$ is a Bernoulli random variable where $Z_{i}=1$ if and only if for the $i$-th sample in $S$ we have $(v\cdot(f(X_{i},p)-p)\geq T)$. Let $\nu_{1}=\nu_{1}(T)=5\exp(-T^{2}/2)/2$, $\nu_{2}=\nu_{2}(T)=\epsilon/(2T^{2})$, and $\nu=\nu(T)=\nu_{1}(T)+\nu_{2}(T)$. We want to prove that

happens with probability at most $\exp(-\Omega(N\epsilon^{2}))$ for any $T>0$.

We have $\mu=\Pr_{Y\sim P}[|v\cdot(f(Y,p)-p)|\geq T]\leq 2\exp(-T^{2}/2)$. Let $T^{\prime}=\Theta(\sqrt{\log(1/\epsilon)})$ be such that $2\mu(T^{\prime})=4\exp(-T^{\prime 2}/2)=\epsilon^{2}/(4T^{\prime 4})=\nu_{2}(T^{\prime})^{2}$. For $T\leq T^{\prime}$, we use bound (i), and for $T\geq T^{\prime}$, we use bound (ii).

Since $\nu\geq\nu_{1}\geq 5\mu/4$, by (i), $\sum_{i}Z_{i}/N\geq\nu$ with probability at most $\exp(-N\ln(5/4)(\nu-\mu)/8)\leq\exp(-N\nu/180)$. When $T\leq T^{\prime}$, $\nu\geq\nu_{2}=\epsilon/(2T^{2})=\Omega(\epsilon/\log(1/\epsilon))$, and so $\sum_{i}Z_{i}/N\geq\nu_{1}(T)$ with probability at most $\exp(-\Omega(N\epsilon/\log(1/\epsilon)))$.

When $T\geq T^{\prime}$, we have $2\mu(T)\leq\nu_{2}(T)^{2}$ and so $\ln(\nu_{2}/\mu)\geq\ln(2/\mu)/2\geq T^{2}/4$. Thus, we have

Using bound (ii), we get that $\Pr[|v\cdot(f(X,p)-p)|\geq T]\geq\nu(T)$ with probability at most $\exp(-\Omega(N\epsilon))$.

In either case, we can take a union bound with the probability that the variance was far above and get that both requirements hold with probability at least $1-\exp(-\Omega(N\epsilon^{2}))$. ∎

For (iii), we will use Claim 25 (ii) with $\epsilon^{\prime}=\epsilon/\ln(d)$. That is, when $N\geq\Omega((m\log(d/\epsilon)+\log(1/\tau))/\epsilon^{\prime 2})=\Omega((m\log(d/\epsilon)+\log(1/\tau))\cdot\log^{2}d\cdot\epsilon^{-2})$, for every $v^{\prime}\in\mathcal{C}$, $T^{\prime}\in\mathcal{T}$, we have

assuming (i) and (ii). This completes the proof of (v).

By a union bound, (i)-(v) all hold simultaneously with probability at least $1-\tau$. ∎

B.2 Omitted Proofs from Section 3.2: Setup and Structural Lemmas

In this section, we prove some structural lemmas that we will need to prove Proposition 9.

First, we note that since the probabilities of the parental configurations are probabilities, the noise will not move them much. Abusing notation, we use $\alpha$ for the empirical minimum parental configuration.

For all $k$, $|\Pr_{S^{\prime}}[\Pi_{k}]-\Pr_{S}[\Pi_{k}]|\leq 2\epsilon$ and $\alpha\geq(C^{\prime}-3)\epsilon\geq\epsilon$.

Proposition 9 requires that $|E|+|L|\leq 2\epsilon|S^{\prime}|$. We have

Since $S$ is $\epsilon$-good, by Lemma 24 (i), $|\Pr_{P}[\Pi_{k}]-\Pr_{S}[\Pi_{k}]|\leq\epsilon$. Since we assume that $\min_{k}\Pr_{P}[\Pi_{k}]\geq 4\epsilon$, $\alpha=\min_{k}\Pr_{S^{\prime}}[\Pi_{k}]\geq\epsilon$. ∎

Our next step is to analyze the spectrum of $M$, and in particular show that $M$ is close in spectral norm to $w_{E}M_{E}$. To do this, we begin by showing that the spectral norm of $M_{S,0}$ is relatively small. Since $S$ is good, we have bounds on the second moments $F(X,p)$. We just need to deal with the error from replacing $p$ with $q$.

Lemma 10. $\|M_{S,0}\|_{2}\leq O(\epsilon+\sqrt{\sum_{k}\Pr_{S}[\Pi_{k}](p_{k}-q_{k})^{2}}+\sum_{k}\Pr_{S}[\Pi_{k}](p_{k}-q_{k})^{2})$.

Let $A_{S}$ denote the second-moment matrix of $(F(X,p)-p)$ under $S$.

First we will show that $M_{S}$ is close to $A_{S}$, and then we will show that their diagonals are close which implies that $M_{S,0}$ is close to $A_{S,0}$.

Now if $|v^{T}(M_{S}-A_{S})v|\leq 4+O(\epsilon)$, then $|v^{T}(M_{S}-A_{S})v|\leq(4+O(\epsilon))\sqrt{\|B\|_{2}}$ and if $||v^{T}(M_{S}-A_{S})v||\geq 4+O(\epsilon)$, then $|v^{T}(M_{S}-A_{S})v|\leq 2\sqrt{|v^{T}(M_{S}-A_{S})v|}\sqrt{\|B\|_{2}}$ and so $|v^{T}(M_{S}-A_{S})v|\leq 4\|B\|_{2}$. Either way, we have $|v^{T}(M_{S}-A_{S})v|\leq O(\max\{\sqrt{\|B\|_{2}},\|B\|_{2}\})$. This holds for all $v$ and so

Now we can bound the spectral norm of $B$ in terms of its Frobenius norm:

Combining this with the bound on $\|M_{S}-A_{S}\|_{2}$ above, we obtain

For the diagonal entries of $M_{S}$ and $A_{S}$, we have

Finally, we can put all this together, obtaining

Lemma 11. $w_{L}\|M_{L}\|_{2}\leq O(\epsilon\log(1/\epsilon)+\epsilon\|p-q\|_{2}^{2})$.

Since $L\subset S$, for any event $A$, we have that $|L|\Pr_{L}[A]\leq|S|\Pr_{S}[A]$. Note that for any $x$, since $((F(x,q)-q)-(F(x,p)-p))_{i}$ is either or $p_{i}-q_{i}$ for any $i$, thus $\|(F(X,q)-q)-(F(X,p)-p)\|_{2}\leq\|p-q\|_{2}$. Since $S$ is $\epsilon$-good for $P$, by Lemma 24 (iii), we have

Also not that $\Pr_{X\in_{u}L}[|v\cdot(F(X,q)-q)|>\sqrt{d}]=0$ since $\|F(X,q)-q\|_{2}\leq\sqrt{d}$. By definition, $\|M_{L}\|_{2}$ is the maximum over unit vectors $v$ of $v^{T}M_{L}v$. For any unit vector $v$, we have We write $f(x)\ll g(x)$ for $f(x)=O(g(x)).$

The last inequality uses $|L|\leq 2\epsilon|S^{\prime}|$ and $|S|\leq(1+2\epsilon)|S^{\prime}|$. ∎

Finally, combining the above results, since $M_{S}$ and $M_{L}$ have small contribution to the spectral norm of $M$ when $\|p-q\|_{2}$ is small, most of it must come from $M_{E}$.

Lemma 12. $\|M-w_{E}M_{E}\|_{2}\leq O(\epsilon\log(1/\epsilon)+\sqrt{\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2}}+\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2})$.

Note that $|S^{\prime}|M=|S|M_{S,0}+|E|M_{E,0}-|L|M_{L,0}$.

Note that each entry of any of these matrices has absolute value at most one since $|F(x,q)-q|_{k}\leq 1$ for all $x\in\{0,1\}^{d}$ and $k\in[m]$. Thus we have

By the triangle inequality, Lemmas 10 and 11, and the assumption that $|E|+|L|\leq 2\epsilon|S^{\prime}|$,

Using Lemma 27, we obtain that $\epsilon\|p-q\|_{2}^{2}\leq\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2}$ and $\sum_{k}\Pr_{S}[\Pi_{k}](p_{k}-q_{k})^{2}\leq\sum_{k}(\Pr_{S^{\prime}}[\Pi_{k}]+2\epsilon)(p_{k}-q_{k})^{2}=O(\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})^{2})$. ∎

B.3 Omitted Proofs from Section 3.2: The Case of Small Spectral Norm

In this section, we will prove that if $\|M\|_{2}=O(\epsilon\log(1/\epsilon)/\alpha)$, then we can output the empirical conditional means $q$.

Since both terms are positive semidefinite, we have

Applying this with $y=q$ and $Y=F(X,q)$ for $X\in_{u}L$ (or $X\in_{u}E$) completes the proof. ∎

When $\Pi_{k}$ does not occur, $F(X,q)_{k}=q_{k}$. Thus, we can write:

By Lemma 27, $|\Pr_{S^{\prime}}[\Pi_{k}]-\Pr_{S}[\Pi_{k}]|\leq 2\epsilon$, so we have $|\Pr_{S^{\prime}}[\Pi_{k}](p_{k}-q_{k})-\Pr_{S}[\Pi_{k}](p_{k}-q_{k})|\leq 2\epsilon|p_{k}-q_{k}|$, i.e., $|z_{k}-z^{\prime}_{k}|\leq 2\epsilon|p_{k}-q_{k}|$, and thus $\|z-z^{\prime}\|_{2}\leq 2\epsilon\|p-q\|_{2}$. ∎

Lemma 15. $\sqrt{\sum_{k}\Pr_{S^{\prime}}[\Pi_{k}]^{2}(p_{k}-q_{k})^{2}}\leq 2\sqrt{\epsilon\|M\|_{2}}+O(\epsilon\sqrt{\log(1/\epsilon)+1/\alpha})$.

Note that $\mu^{S^{\prime}}=q$. By Lemma 14, $\|(\mu^{S}-q)-(D^{2}(p-q))\|_{2}\leq O(\epsilon(1+\|p-q\|_{2}))$. Recall that $w_{L}=|L|/|S^{\prime}|$ and $w_{E}=|E|/|S^{\prime}|$. By the triangle inequality,

where we used the assumption that $\alpha/\epsilon$ is at least a sufficiently large constant and that $\|D^{-1}\|_{2}=1/\sqrt{\alpha}$. When this last term is smaller than $\|D^{2}(p-q)\|_{2}/8$, rearranging the inequality gives that

Otherwise, we have $\|D^{2}(p-q)\|_{2}=O\left(\sqrt{\epsilon\|D^{2}(p-q)\|_{2}/\sqrt{\alpha}}\right)$, and so $\|D^{2}(p-q)\|_{2}=O(\epsilon/\sqrt{\alpha})$. In either case, we obtain

Corollary 16 (Part (i) of Proposition 9). If $\|M\|_{2}\leq O(\epsilon\log(1/\epsilon)/\alpha)$, then

Recall that $\alpha=\min_{k}|\Pr_{S^{\prime}}[\Pi_{k}]$. By Lemma 27, $|\Pr_{S^{\prime}}[\Pi_{k}]-\Pr_{S}[\Pi_{k}]|\leq 2\epsilon$ for any $k$. Since $S$ is $\epsilon$-good, $|\Pr_{P}[\Pi_{k}]-\Pr_{S}[\Pi_{k}]|\leq\epsilon$. Combining these, we obtain $|\alpha-\min_{k}\Pr_{P}[\Pi_{k}]|\leq 3\epsilon$. By assumption $\min_{k}\Pr_{P}[\Pi_{k}]\geq 4\epsilon$ , so we have $\alpha=\Theta(\min_{k}\Pr_{P}[\Pi_{k}])$.