Year: 2019

Abstract:
Feature selection is a pervasive problem. The discovery of relevant features can be as important for performing a particular task (such as to avoid overfitting in prediction) as it can be for understanding the underlying processes governing the true label (such as discovering relevant genetic factors for a disease). Machine learning driven feature selection can enable discovery from large, high-dimensional, nonlinear observational datasets by creating a subset of features for experts to focus on. In order to use expert time most efficiently, we need a principled methodology capable of controlling the False Discovery Rate. In this work, we build on the promising Knockoff framework by developing a flexible knockoff generation model. We adapt the Generative Adversarial Networks framework to allow us to generate knockoffs with no assumptions on the feature distribution. Our model consists of 4 networks, a generator, a discriminator, a stability network and a power network. We demonstrate the capability of our model to perform feature selection, showing that it performs as well as the originally proposed knockoff generation model in the Gaussian setting and that it outperforms the original model in non-Gaussian settings, including on a real-world dataset.

heading: 1 INTRODUCTION
text: Feature selection is a pervasive problem. Often the goal is to discover features that are relevant to a particular outcome, either for the sake of discovery itself or to aid in prediction [16; 25]. When the focus is on discovery, feature selection methods typically focus on trying to control either the Family-Wise Error Rate (FWER) or the False Discovery Rate (FDR). The FWER measures the probability of making a single false discovery (a Type I error) among the selected features (i.e. selecting one which is not relevant), whereas the FDR measures the proportion of false discoveries made (i.e. the proportion of selected features which are false). Controlling FWER, however, leads to reduced power (i.e. selecting fewer relevant variables) since it controls the probability of making any false discovery, whereas FDR tries to control the proportion of false discoveries.
Controlling the FDR is important [5; 6; 3]. Often, data-driven feature selection will be used to select a set of candidate features for further investigation. When further investigation is expensive (for example when further investigation would involve conducting new experiments and collecting more data), a method that cannot control the FDR may result in a large amount of wasted resources, with no guarantee that anything meaningful will be discovered. On the other hand, being able to control the FDR at, say, 10% ensures that at most, 10% of the spent resources are wasted, and 90% are in fact spent on discovering positive, useful results. It should be noted, however, that estimating
the FDR of a method empirically is hard in practice, since we do not have access to the ground truth relevance. As such, a theoretical analysis of the method and its (potential) FDR-controlling properties must be carried out, which does not exist for many existing feature selection methods.
[3] is the seminal paper on the knockoff framework, which is an innovative FDR-controlling feature selection method. Knockoffs are features that are generated to “look like” the real features but be conditionally independent of the label given the real features. Feature statistics (such as the coefficients of a LASSO [32]) are compared between the real features and their knockoffs and a selection is made when this difference is sufficiently large. Performing the selection in this way allows for an estimate of the FDR to be obtained and the selection threshold can be adjusted to control the FDR at the selected level. In the original paper, the relationship between the label and the features is constrained to be of a very specific form; in [7], they remove this constraint and instead provide a theoretical analysis that shifts the burden of knowledge onto knowing the underlying feature distribution. Unfortunately, while the theoretical results hold for any feature distribution, they rely on being able to generate valid knockoffs, for which [7] only provide a method for generating knockoffs when the distribution is a (known) multivariate Gaussian distribution. In this paper, we modify the Generative Adversarial Networks (GAN) [11] framework to address this problem, allowing us to generate knockoffs for any distribution (and without any prior knowledge of it). GANs have been shown to be a powerful method for learning to generate complex distributions [24; 20; 2].
Our main contribution is in modifying the discriminator used in the GAN framework in such a way that the generator learns to generate knockoffs satisfying the necessary swap condition [7] which requires that when a feature and its knockoff are swapped, the joint distribution remains unchanged. In addition, we propose a method for maximizing the power of our model using Mutual Information Neural Estimation (MINE) [4] and investigate a regularization method to improve the stability of training. Our model consists of four networks: (1) a generator network that takes as input noise and the real features, and outputs a set of candidate knockoff features; (2) a discriminator network taking as input “swapped” feature-knockoff features that attempts to determine which variables have been swapped; (3) a Wasserstein GAN discriminator that we use as a regularization term; and (4) a MINE network that estimates the mutual information between each feature-knockoff pair allowing us to maximize the power of the knockoff procedure.

heading: 2 RELATED WORKS
text: Feature selection is a well-studied problem with a wealth of related works ([12; 31; 17; 22] provide a summary of a lot of existing literature); however, most methods do not attempt to control the FDR. The most common feature selection method for FDR control is the Benjamini-Hochberg (BHq) procedure and its variants [5; 6], which relies on obtaining valid marginal p-values for each selection.
Knockoffs are an active area of research [9; 18; 10]. The notion of a knockoff was first introduced in [3] with the theory there requiring that the relationship between the features and the label be of a specific form. In [7], they build on the knockoff framework, removing this requirement but instead shifting the requirement to one of knowing the distribution of the features. As noted in the introduction, the theory in [7] holds independent of the distribution of the features - relying only on being able to generate valid knockoffs (which exist for any distribution of features). However, they only propose a method for generating knockoffs when the distribution of features is jointly Gaussian. While they do propose a method for generating approximate knockoffs in the non-Gaussian setting (by simply approximating the features as Gaussian), the guarantees on FDR control do not hold for their approximate knockoffs. In [26] and [10], they add to the class of constructible knockoffs, describing methods for constructing knockoffs for Markov Chains, Hidden Markov Models and Gaussian Mixture Models. Though once again, knowledge of the full distribution is still necessary for their construction.
In this paper we use a framework motivated by GANs [11] to learn to generate knockoffs without any assumptions on the distribution of the features. To do this, we modify the discriminator so that rather than trying to determine whether a sample is real or fake, it attempts to identify which components have been “swapped”. In [38], an unconventional discriminator is used that performs component-wise discrimination for the purpose of imputation. While the problem addressed in that paper is different to the one here, the key idea relies on a similar modification to the discriminator to be able to appropriately guide the generator.
In order to maximize the power of our variable selection mechanism, it will be desirable that the feature-knockoff pairs are ”as independent as possible” (this is discussed in [7]). In order to achieve this we will investigate the use of a promising recent paper, MINE [4]. MINE proposes a neural architecture and training procedure capable of estimating the mutual information between two random variables. As the mutual information between two random variables is zero only when they are independent, we will use this as a measure of independence and attempt to minimize it during the training of our modified GAN.

heading: 3 BACKGROUND
text: In this section we introduce our notation and define knockoffs as in [7]. Let us denote the feature space by X and the label space by Y . Let the dimension of X be d. Suppose that X = (X1, ..., Xd) and Y are random variables over X and Y . As in [7], we will work with the notion of a null set. Definition 1. A variable Xj is said to be ”null” if and only if Y is independent of Xj conditional on {Xi : i 6= j}. We defineH0 to be the set of all null variables.
Our goal will be to discover as many relevant features as possible while controlling the FDR. For a given (potentially random) selection procedure that selects Ŝ ⊂ {1, ..., d}, we define the FDR to be
FDR = E [ |Ŝ ∩ H0| |Ŝ| ] .
Note that this agrees with the usual notion of FDR (i.e. when defined in terms of the Markov blanket) under mild assumptions (for a more thorough discussion see [7]).

heading: 3.1 KNOCKOFFS
text: Definition 2. A knockoff [7] for X is a random variable X̃ ∈ X satisfying the following two properties:
(X, X̃) d. = (X, X̃)swap(S) (1)
X̃ ⊥⊥ Y |X (2) for all S ⊂ {1, ..., d} where (·, ·)swap(S) denotes the vector obtained by swapping the ith component with the (i+ d)th component for each i ∈ S and d.= is equality in distribution.
In order to use knockoffs for feature selection, we must define an appropriate feature statistic, Wj , that depends on X, X̃ and Y , i.e. Wj = wj((X, X̃), Y ) for some function wj . This function wj must satisfy the following flip-sign property:
wj((X, X̃)swap(S), Y ) = { wj((X, X̃), Y ) if j /∈ S −wj((X, X̃), Y ) if j ∈ S.
(3)
One of the procedures used in [7] to construct these statistics is to perform LASSO, treating the augmented feature-knockoffs as the features on which to regress. This gives LASSO coefficients b1, ..., b2d, and the statistic Wj is set to be the LASSO Coefficient Difference given by
Wj = |bj | − |bj+d|.
Note that the FDR control guarantees hold independently of the choice of statistic, but a poorly chosen statistic can significantly impact the power of the test. In particular, using the LASSO Coefficient Difference in non-linear settings can yield few discoveries. The focus of this paper, however, is on generating the knockoffs, not on the statistic used on top of the generated knockoffs and so in our synthetic experiments, we use a linear model for Y to be able to draw fair comparisons between our model and [7]. In the real world data experiment, we use a statistic based on Random Forests for both methods [37].
The following result from [7] depends only on having obtained knockoffs that satisfy Definition 2 and feature statistics satisfying (3) (and in particular do not depend specifically on using LASSO to obtain the statistics).
Theorem 1. Let q ∈ [0, 1]. Given test statistics, W1, ...,Wd, satisfying (3), let
τ = min { t > 0 :
1 + |{j : Wj ≤ −t}| |{j : Wj ≥ t}|
≤ q } .
Then the procedure selecting the variables
Ŝ = {j : Wj ≥ τ} controls the FDR at level q, i.e.
E [ |Ŝ ∩ H0| |Ŝ| ∨ 1 ] ≤ q.

heading: 4 KNOCKOFFGAN
text: It should be noted that in order to satisfy equation (2), it simply needs to be the case that the knockoffs are constructed without looking at the label, Y . In order to satisfy equation (1) we use a modified GAN framework, which gives us the flexibility to learn to generate knockoffs without any assumptions on the distribution of the original features.

heading: 4.1 GENERATOR
text: The generator, G, will be a function G(·, ·;φ) : X × [0, 1]c → X , parametrized by φ that takes a realization x of X and random noise, z ∼ U([0, 1]c), as inputs and outputs knockoff features x̃. We define X̃ := G(X, z). We model G as a fully connected neural network with weights φ.

heading: 4.2 DISCRIMINATOR
text: The main innovation of our paper is in defining the discriminator. Equation (1) imposes a condition on the joint distribution of (X, X̃) and as such we must define a discriminator with a loss that is (not necessarily uniquely) minimized only for joint distributions satisfying this condition. To that end, the discriminator, D, will be a function D(·;ψ) : X × X → [0, 1]d that takes as input a swapped sample-knockoff pair (x, x̃)swap(S) and outputs a vector in [0, 1]d with the ith component of D((x, x̃)swap(S)) corresponding to the probability that i ∈ S. The discriminator is attempting to detect which variables have been swapped and, intuitively, when the discriminator is unable to determine this, the swapped and unswapped joint distributions must be the same.
The loss we use to train the discriminator is the multi-output cross-entropy loss given by LD = ∑
S∈{0,1}d EX∼PX [EX̃∼P̃X(X)[S·log(D((X, X̃)swap(S)))+(1−S)·log(1−D((X, X̃)swap(S)))]]
(4)
where · is the standard dot, 1 = (1, ..., 1), S = (S1, ..., Sd) with Si = I(i ∈ S) (I is the indicator function) and log is taken element-wise. The following theorem is our main theoretical result, which states that the training regime employed by KnockoffGAN will result in a procedure that generates valid knockoffs. Theorem 2. Equation (4) is maximized (with respect to G) if and only if equation (1) is satisfied by G.
Proof. The proof, alongside supporting theoretical results, can be found in the Appendix.
In practice, the sum is too computationally expensive (O(2d)) to calculate and so we perform stochastic gradient descent using minibatches with S sampled uniformly from {0, 1}d, independently for each sample in the minibatch.
We also found that training with respect to the full loss resulted in a poor performance, particularly when d is large. We found that the discriminator struggled to learn anything when asked to find the full swap vector, and the poor discriminator resulted in a poorly trained generator. In order to overcome this, we introduce a hint vector - first introduced in [38] - that we use to reveal partial information to the discriminator about the swap vector. We do this by using the hint to reveal some, but not all, of the components of S to the discriminator. In doing so, we reduce the burden of the discriminator from needing to determine the entire swap vector to only needing to determine some of the swap vector.
Formally, the hint, H, will be a random variable depending on S, that we pass to the discriminator, alongside (X, X̃)swap(S). We use the hint to control the amount of information we pass to D about S before asking D to predict S. In practice, our hinting mechanism involves sampling a multivariate Bernoulli random variable, B from i.i.d. components, which each take value 1 with probability 0.9. The hint is then constructed by setting Hi = Si if Bi = 1 and Hi = 0.5 if Bi = 0. The discriminator is therefore being asked only to predict the values of S for which Bi = 0; the others, D is able to directly infer from Hi. In order to avoid overfitting to the hint, it becomes necessary to remove these terms from our loss. Our loss now becomes
LD = ∑
S∈{0,1}d EX∼PX
[ EX̃∼P̃X(X) [ EH∼PH|S [ (S (1−B)) · log(D((X, X̃)swap(S),H)) (5)
+ ((1− S) (1−B)) · log(1−D((X, X̃)swap(S),H)) ]]]
where denotes element-wise multiplication and the expectation over B is implicit in the expectation over H.

heading: 4.3 STABILITY
text: We found that adding a regularization term in the form of a Wasserstein GAN discriminator (with Gaussian Process (GP) regularization) [2], f , aided performance. We note that when equation (1) holds, we must have that X d.= X̃ and so the addition of this regularizing term does not affect the optimal solution to our loss. We model f as a fully connect neural network with weights ν. The loss is given by
Lf = E [ f(X)− f(X̃)− η(||∇X̂f(X̂)||2 − 1) 2 ]
where ∼ U [0, 1], X̂ = X + (1 − )X̃ and η is a hyper-parameter (set to 10 in practice). Note that we have rewritten the loss to be the negative of the one given in [2], allowing us to write our overall objective as a minimax problem. This loss is added to the generator loss as an additional regularization term.

heading: 4.4 MAXIMIZING POWER
text: As noted in [7], it is intuitive that in order to maximize the power of the knockoff selection procedure, we wish to make Xj and X̃j as ”independent” as possible. Doing so ensures that as little as possible of the dependence between the real feature and the label is present between the knockoff and the label; this allows us to determine whether or not the relationship between the feature and label is only through the feature’s correlation with other features, or is in fact a true signal.
In order to achieve maximal independence, we look to minimize the mutual information between each feature and its knockoff. Actually computing the true mutual information requires access to both the joint density of the feature-knockoff pairs and to the marginal densities of each feature and knockoff, which we do not have.
Instead, we look to a promising recent work, Mutual Information Neural Estimation (MINE [4]), that provides a framework for estimating the mutual information using neural networks. To do so, they estimate the mutual information between random variables U and V by performing gradient ascent on the following objective:
sup θ∈Θ EP(n)UV [Tθ]− log(EP(n)U ⊗P(n)V [e Tθ ])
where PUV denotes the joint measure of (U, V ) with PU = ∫ V dPUV and PV = ∫ U dPUV denoting the marginal measures. (n) denotes the empircal distribution associated with n i.i.d. samples.
Using MINE we approximate the mutual information between each pairXj and X̃j by using d neural networks1, T 1, ..., T d, each parametrized by θ1, ..., θd, that we refer to collectively as the power network, and will write P to denote the collection of networks T 1, ..., T d. The mutual information is added using a trade-off parameter λ to the loss for G. Formally, define LP by
LP = d∑ j=1 ( n∑ i=1 (T jθj (x (i) j , x̃ (i) j ))− log( n∑ i=1 exp(T jθj (x (κ(i)) j , x̃ (i) j ))) )
where κ is a random permutation of [n]2 and (i) denotes the ith sample - noting that dependence on G is through X̃.

heading: 4.5 FINAL OBJECTIVE
text: The resulting minimax game played by G, D, W and P is given by
min G ( max D (LD) + λmax P (LP ) + µmax f (Lf ) )
where λ, µ are hyper-parameters (set to 1 in the experiments section).
We train each of G, D, W and P iteratively. Pseudo-code of our knockoff construction algorithm can be found in Algorithm 1 and a visual representation of our architecture in Fig. 1.
After generating knockoffs, feature statistics are computed according to some procedure (in the synthetic experiments we use LASSO and in the real data experiment we use a Random Forestbased statistic [37]). Features are then selected based on these statistics according to Theorem 1.

heading: 5 EXPERIMENTS
text: In this section we demonstrate the capability of our method to match the results of [7] in settings where their model is correctly specified (i.e. when the underlying feature distribution is Gaussian) and then go on to show, in settings where the underlying feature distribution is non-Gaussian, that our method is able to outperform their Gaussian approximation. We compare to two versions of the BHq method [5; 6] to provide a baseline.
We also perform a qualitative analysis of KnockoffGAN on a real-world dataset. We compare features found by KnockoffGAN to PubMed literature and show that KnockoffGAN discovers several meaningful features for 2 different disease outcomes.
1In practice we use a single neural network with diagonalized weights to parallelize these networks. 2In our pseudo-code, we write U(Sn) to denote the uniform distribution over the set of all permutations of
[n] = {1, ..., n}.
Algorithm 1 Pseudo-code of KnockoffGAN 1: Inputs: mini-batch size nmb > 0, Initialize parameters φ, ψ, ν, θ1, ..., θd 2: while Converge do 3: Discriminator Update 4: Sample x1, ...,xnmb from D, z1, ..., znmb ∼ PZ 5: Sample S1, ...,Snmb
i.i.d∼ U({0, 1}d), b1, ...,bnmb ∼ Ber(0.9) 6: for i = 1, ..., nmb do 7: x̃i ← G(xi, zi;φ) 8: hi = Si bi + 0.5(1− bi) 9: Update D by ascending its stochastic gradient
∇ψ nmb∑ i=1 [ (Si (1− bi)) · log(D((xi, x̃i)swap(S)),hi)
+ ((1− Si) (1− bi)) · log(1−D((xi, x̃i)swap(S),hi)) ]
10: MINE Update 11: Sample x1, ...,xnmb from D, z1, ..., znmb ∼ PZ , κ ∼ U(Snmb) 12: for i = 1, ..., nmb do 13: x̃i ← G(xi, zi;φ) 14: for j = 1, ..., d do 15: Update Tj by ascending its stochastic gradient
∇θj (∑nmb i=1 T j θj (x (i) j , x̃ (i) j ) ) − log (∑nmb i=1 exp(T j θj (x (i) j , x̃ (κ(i)) j )) ) 16: WGAN-GP Update 17: Sample x1, ...,xnmb from D, z1, ..., znmb ∼ PZ 18: for i = 1, ..., nmb do 19: Sample ∼ U [0, 1] 20: x̃i ← G(xi, zi;φ) 21: x̂i = xi + (1− )x̃i 22: Update f by ascending its stochastic gradient
∇ν ∑nmb i=1 [ f(xi)− f(x̃i)− η(||∇x̂if(x̂i)||2 − 1)2 ] 23: Generator Update 24: Sample x1, ...,xnmb from D, z1, ..., znmb ∼ PZ 25: Sample S1, ...,Snmb
i.i.d∼ U({0, 1}d), κ ∼ U(Snmb) 26: for i = 1, ..., nmb do 27: x̃i ← G(xi, zi;φ) 28: Update G by descending its stochastic gradient
∇φ(LD + λLP + µLf )

heading: 5.1 SYNTHETIC DATA EXPERIMENTS

heading: 5.1.1 SIMULATION SETTINGS
text: Evaluating feature selection methods on real data is difficult as we do not have access to the ground truth. To evaluate KnockoffGAN, we conduct a series of experiments using synthetic data, replicating those carried out in [7] and extending them to more general settings. In each of the following synthetic experiments, we set the feature dimension to be d = 1000 and the number of samples to be n = 3000. For each feature distribution we perform two experiments:
1. Y-Logit: P (Y = 1|X) = exp(m(X))(1+exp(m(X)))
2. Y-Gaussian: Y ∼ N (m(X), 1)
where m(X) = ∑60 i=1 αδiXi with δi ∈ {−1, 1} sampled uniformly and then fixed for each experiment. α controls the strength of the influence that X has on Y , and in the experiments we vary this (as in [7]). Note that for the auto-regressive settings (found in Section 5.1.2 and the Appendix) the relevant variables are sampled uniformly at random from among the 1000 features (rather than being the first 60); in the non-auto-regressive settings this is not necessary.
We report the True Positive Rate (TPR), which is also commonly referred to as the power, defined as
TPR = |Ŝ ∩ S∗| |S∗|
(6)
where S∗ = {1, ..., d} \ H0 is the set of all non-null features. We also report the FDR to verify that the methods do indeed control it at the specified level which we set to be 10%. Note that we are not using FDR as a metric - a lower FDR is not desirable when we set the threshold to 10%. In fact, we want the methods to be as close to 10% as possible (so that they are achieving maximum power). We perform 100 replications of each experiment and report the average TPR and FDR.

heading: 5.1.2 GAUSSIAN SETTINGS
text: We begin by replicating the setup from [7] in which the underlying feature distribution is Gaussian. In this setting, we do not expect KnockoffGAN to perform better than the original knockoff framework as the original framework assumes a Gaussian distribution. Our goal here is simply to achieve a similar performance, demonstrating that little performance is lost even when the distribution is known to be Gaussian.
In the first experiment that we replicate from [7], the features are set to be auto-regressive (AR(1)) with Gaussian marginal distributions, i.e. Xi = φXi−1 +Zi with Zi being chosen such thatXi
i.i.d.∼ N (0, 1n ). In this experiment we vary φ, which determines the correlation between features, rather than α. We fix α = 3.5 for Y-Gaussian and α = 10 for Y-Logit. The results are reported in Fig. 2.
As in [7], we observe that BHq Marginal, which tests for marginal independence of the feature from Y , suffers from severely increased FDR as we increase the correlation, invalidating the seemingly good TPR. To make the remaining results clearer, we omit BHq marginal from the rest of this section. Aside from this, we see in Fig. 2, that the other methods control the FDR at or very close to the specified 10% threshold. We also see that across the entire range of α, KnockoffGAN achieves a very similar TPR to the original Knockoff framework.
In the second experiment, we set the underlying feature distribution to be i.i.d. Gaussian. We found in this case also that KnockoffGAN was able to control the FDR and achieve a similar TPR to the
original knockoff framework. More details of this experiment and the results for it can be found in the Appendix

heading: 5.1.3 NON-GAUSSIAN SETTINGS
text: We now move on to the key results for the paper in which the underlying feature distribution is no longer Gaussian. In this setting, we expect to outperform the original Knockoff framework due to the fact that they approximate the distribution as Gaussian. In particular, when this approximation is poor, the knockoffs are no longer valid and as such no FDR guarantees can be given. On the other hand, KnockoffGAN does not place any requirements on the distribution of the features and as such is able to generate valid knockoffs.
We performed experiments for several different underlying feature distributions, and found that KnockoffGAN achieved a higher TPR than the original knockoff framework in all cases, while controlling the FDR at the specified level. We give the results for X coming from a 4-Gaussian mixture model in Fig. 3 - results for Uniform, Dirichlet, and other (2 and 3) Gaussian mixture models can be found in the Appendix.
To create our 4-mixture model, we set the means (m1,m2,m3,m4) of the 4 Gaussians to be:
• m1i = 1 for i = 1 to 100 and 0 for i = 101 to 1000, • m2i = 1 for i = 1 to 50 and -1 for i = 51 to 100 and 0 for i = 101 to 1000, • m3i = −1 for i = 1 to 50 and 1 for i = 51 to 100 and 0 for i = 101 to 1000, • m4i = −1 for i = 1 to 100 and 0 for i = 101 to 1000.
We scale the variance of each Gaussian to be such that the overall variance of each feature is 1n .
We see in Fig. 3 that KnockoffGAN consistently outperforms the original knockoff framework, achieving a higher TPR across the entire range of α while consistently controlling the FDR at 10%. In fact, in the Y -Gaussian setting we see that the original knockoff framework performs almost identically to BHq Maximum Likelihood.

heading: 5.1.4 IMPACT OF WGAN REGULARIZATION
text: We conclude the synthetic experiments by demonstrating the effect of the WGAN regularizer3. We conduct this experiment using an auto-regressive model with U(− √ 3/n, √ 3/n) marginal distributions. We fix α = 5 for Y-Logit and α = 2.5 for Y-Gauss.
3The WGAN regularizer was included in all previous experiments.
As we see in Fig. 4, the WGAN regularizer has a significant effect on the results, with the improvement in some places being almost as much as KnockoffGAN without WGAN makes over the original knockoff framework. As noted in Section 4.3, there is no trade-off introduced by the inclusion of this regularizer; the optimal solution to the loss is unchanged and therefore this regularization is ”free” in terms of FDR control, but as demonstrated improves TPR performance.

heading: 5.2 REAL DATA EXPERIMENT
text: In this section we use a biobank dataset (with 387 dimensions) to qualitatively analyze the performance of KnockoffGAN. We use KnockoffGAN to select features for two different outcomes: (1) Cardiovascular Disease (CVD) and (2) Diabetes and then use PubMed literature to asses the validity of the selected features.
We found that the original knockoff framework was unable to select even the most well-known features (such as Age and Sex for CVD [15]), even when the FDR threshold was increased to 20%. Therefore, there are no relevant features to report for the original knockoff framework and so Table 1 contains only the features selected by KnockoffGAN that were deemed relevant by PubMed literature. For this the FDR threshold was set to 5% so that the number of discoveries was manageable for cross-reference with PubMed.
As we see in Table 1, KnockoffGAN discovers 9 relevant features for CVD and 8 relevant features for diabetes. Some of the relevant features, such as Age, Sex and Long-standing illness for CVD are
well-known relevant features. The remaining features are all supported by the literature in PubMed. While this is a qualitative result (it relies on using PubMed as the ground truth), we do believe this demonstrates that KnockoffGAN is a significant improvement over the original knockoff generation procedure.

heading: 6 CONCLUSION
text: In this paper we built on the knockoff framework introduced in [3] by developing a novel GAN framework, KnockoffGAN, capable of generating knockoffs with no assumptions on the underlying data. We demonstrated through a series of experiments on a range of synthetic datasets and on a real world dataset that our method improves on the performance of the original knockoff framework.
While we feel this is a significant step towards being able to generate knockoffs for any data, there is still more work to be done. In particular, generalizing this method to time-series data would be non-trivial, and would be an interesting avenue for further investigation.

heading: ACKNOWLEDGEMENT
text: The authors would like to thank the reviewers for their helpful comments. The research presented in this paper was supported by the Office of Naval Research (ONR) and the NSF (Grant number: ECCS1462245, ECCS1533983, and ECCS1407712).

heading: APPENDIX

heading: THEORETICAL RESULTS
text: In order to prove Theorem 2, we use similar techniques to those used in the original GAN paper [11]. In what follows, we analyze the minimax game defined by:
min G max D ∑ S∈{0,1}d EX∼PX [EX̃∼P̃X(X)[S·log(D((X, X̃)swap(S)))+(1−S)·log(1−D((X, X̃)swap(S)))]] (7) where we have used the version of LD given by equation (4) in the main manuscript (i.e. without hinting). The theoretical results that follow are proven only for this version of the loss, though we do believe that the theorem holds more generally for the hinting version of the loss - this is backed up by our empirical results demonstrating strict FDR control while using the hint mechanism. After proving that the optimal solution to this game does indeed provide us with valid knockoffs, we will then show that the additional term Lf does not change the optimal solution. Let p be the density of (X, X̃).
We begin by stating a lemma, that follows from a similar proof to Proposition 1 in [11].
Lemma 1. Let (x, x̃) ∈ X × X . Then for a fixed generator, G, the ith component of the optimal discriminator, D∗((x, x̃)) is given by
D∗((x, x̃))i = p((x, x̃)swap({i}))
p((x, x̃)swap({i})) + p((x, x̃)) (8)
for each i ∈ {1, ..., d}.
Proof. The proof of this involves some basic integral manipulation to get that the objective in (7) can be rewritten as
d∑ i=1 ∫ X×X logD((x, x̃))ip((x, x̃)swap({i})) + log(1−D((x, x̃))i)p((x, x̃))dx.
We then observe that y 7→ a log y + b log(1 − y) achieves its maximum in [0, 1] at aa+b and so the objective is maximized (with respect to D, for fixed G) when
D∗((x, x̃))i = p((x, x̃)swap({i}))
p((x, x̃)swap({i})) + p((x, x̃))
for each i ∈ {1, ..., d}.
With this lemma, we are now able to prove our key result.
Theorem 3. Equation (7) is maximized (with respect to G) if and only if equation (1) (in the main paper) is satisfied by G.
Proof. We begin by rewriting our loss, substituting in D∗, to give us the following loss for G:
LG = ∑
S⊂{1,...,d} EX,X̃ (∑ i∈S log p((x, x̃)swap(S\i)) p((x, x̃)swap(S\i)) + p((x, x̃)swap(S))
+ ∑ i/∈S log p((x, x̃)swap(S)) p((x, x̃)swap(S∪i)) + p((x, x̃)swap(S)) ) where we note that
((x, x̃)swap(S))swap({i}) = { (x, x̃)swap(S\i) if i ∈ S (x, x̃)swap(S∪i) if i /∈ S.
Then by inspecting each term in the sum, we see that each term is a KL-divergence term that is minimized only when p((x, x̃)swap({S\i})) = p((x, x̃)swap(S)) (9) and
p((x, x̃)swap({S∪i})) = p((x, x̃)swap(S))
for every i ∈ {1, ..., d}, every S ⊂ {1, ..., d} and each (x, x̃) ∈ X × X . By iteratively applying equation (9), we get that
p((x, x̃)swap(S)) = p((x, x̃)swap(S\1)) = ... = p((x, x̃)swap(S\{1,...,d−1}))
= p((x, x̃)swap(S\{1,...,d})) = p((x, x̃)
proving the theorem.
Lemma 2. The addition of the term Lf to our loss, does not affect the optimal solution to it.
Proof. By theorem 2, it suffices to show that any distribution satisfying equation (1) also minimizes maxf Lf . But we note that, as shown in [2], maxf Lf (or rather supf Lf ) is the Wasserstein distance between X and X̃, which is 0 (and minimal) when X d.= X̃. It therefore suffices to show that equation (1) implies X d.= X̃.
Let S = {1, ..., d}. Then if (X, X̃) satisfy equation (1), we get that (X, X̃) d.= (X̃,X). Since the joint distributions are equal, it follows that the marginal distributions are equal and so by projecting onto the first d variables, we get that X d.= X̃.

heading: MINE
text: We state the key theory used by MINE to estimate the mutual information. For full details see the original paper, [4].
The mutual information is defined as
I(U ;V ) = ∫ U×V log dPUV dPU ⊗ dPV dPUV
where U and V are random variables over some spaces U and V , respectively with joint measure PUV and marginal measures PU = ∫ V dPUV and PV = ∫ U dPUV , respectively.
The mutual information can also be characterized by the Kullback-Leibler divergence, DKL, as
I(U ;V ) = DKL(PUV ||PU ⊗ PV )
The Donsker-Varadhan representation then gives us for any two probability measures P and Q over a probability space Ω
DKL(P||Q) = sup T :Ω→R EP[T ]− log(EQ[eT ])
where the supremum is taken over all functions T such that the two expectations are finite.
A simple corollary of this is that fixing a classF of functions (such as a parametrized class {Tθ : θ ∈ Θ}) over which the supremum is taken will provide us with a lower bound for the mutual information that approaches the true mutual information as the class becomes sufficiently rich. MINE [4] fix the class to be parametric in this way - given a fixed neural network architecture, they let F = {Tθ : θ ∈ Θ} be the set of all functions parametrized by this network.

heading: IMPLEMENTATION OF KNOCKOFFGAN
text: In the experiments, the depth of the generator, discriminator and WGAN-GP networks is set to 4 and power network is set to 3. The number of hidden nodes in each layer is d/4, d/16 d/4 for the generator, discriminator, and WGAN-GP, respectively. For the power network, we use 2 diagonal matrices for each layer to make two hidden nodes for each feature separately. We use ReLu and tanh as the activation functions for each layer except for the output layer where we use a linear activation function for the generator, power network and WGAN-GP networks and sigmoid activation function for the discriminator network. The number of samples in each mini-batch is 128. KnockoffGAN is implemented in tensorflow.

heading: DETAILS OF BENCHMARKS
text: We use the following links for the implementations of 3 benchmarks.
• Knockoff: http://web.stanford.edu/group/candes/knockoffs/ software/knockoff/index.html
• BHq Max: http://web.stanford.edu/group/candes/knockoffs/ software/knockoff/tutorial-4-r.html
• BHq Marginal: Modifying the original code of BHq Max in http://web.stanford. edu/group/candes/knockoffs/software/knockoff/tutorial-4-r. html
Except for the knockoff generation step, KnockoffGAN follows the same procedures as the original knockoff framework described at http://web.stanford.edu/group/candes/ knockoffs/software/knockoff/tutorial-2-r.html to select features.

heading: ADDITIONAL EXPERIMENTS

heading: INDEPENDENT GAUSSIANS
text: In the following experiment, features were taken to be i.i.d. Gaussian, with mean 0 and variance 1n , i.e. X ∼ N (0, 1nIn). The results are reported in Fig. 5.
As we see in Fig. 5, all methods control the FDR at or very close to the specified 10% threshold. We also see that across the entire range of α, KnockoffGAN achieves a very similar TPR to the original knockoff framework.

heading: NON-GAUSSIAN SETTINGS

heading: INDEPENDENT UNIFORM
text: In this experiment we set the feature distribution to be a Uniform distribution with mean 0 and variance 1n (to be consistent with the Gaussian experiments). Fig. 6 displays the results for each component of X being i.i.d. U(− √ 3/n, √ 3/n). Once again we see that KnockoffGAN consistently outperforms the original knockoff framework, achieving a higher TPR across the entire range of α in both settings.
3/n). TPR is used to
quantify performance and FDR is reported to verify that it is at the specified threshold (10%).

heading: DIRICHLET
text: In this experiment we set the feature distribution to be a Dirichlet(1, ..., 1) distribution - i.e. the uniform distribution over the (d− 1)-simplex. Correlation here exists through the requirement that∑d i=1Xi = 1.

heading: GAUSSIAN MIXTURE MODELS
text: For the GMM2 model we set the means (m1,m2) of the 2 Gaussians to be:
• m1i = 1 for i = 1 to 100 and 0 for i = 101 to 1000, • m2i = −1 for i = 1 to 100 and 0 for i = 101 to 1000.
For the GMM3 model we set the means (m1,m2,m3) of the 3 Gaussians to be:
• m1i = 1 for i = 1 to 100 and 0 for i = 101 to 1000, • m2i = 1 for i = 1 to 50 and -1 for i = 51 to 100 and 0 for i = 101 to 1000, • m3i = −1 for i = 1 to 100 and 0 for i = 101 to 1000.

heading: HYPER-PARAMETER ANALYSIS
text: Hyper-parameter selection in the feature selection problem is difficult; hyper-parameter selection cannot be performed using cross-validation as we do not have access to ground truth. The hyperparameters must therefore be fixed a priori. For this we believe that λ = 1 and µ = 1 are perhaps the most canonical choice we could make and thus these were used in the main manuscript. In the following experiments, we investigate the sensitivity of the results to various settings of λ and µ. The results below are in the auto-regressive Uniform setting with Y-Logit and FDR set to 10%.
As can be seen in Table 2 and 3, the performance of KnockoffGAN is not sensitive to the value of λ and µ. The only significant difference can be seen when either λ or µ is set to 0 which represents the exclusion of either the power network or WGAN-discriminator network, respectively, from the model. In particular, the lack of sensitivity to µ aligns with Lemma 2 in which we see that there is not a trade-off between Lf and LD but rather the two can be simultaneously minimised. To further understand the effects of the WGAN network, we report the final values of the other losses (LD and LP ) when µ = 1 and µ = 0 (i.e. with and without the WGAN regularisation). The results are given below in the auto-regressive Uniform setting with Y-Logit and FDR set to 10%.
As can be seen in Table 4, the inclusion of the WGAN network improves the final value of LD, bringing it significantly closer to its optimal value of log(2) ≈ 0.693. On the other hand, the power network loss is increased, however, this trade-off is expected and acceptable - the FDR guarantees rely on the discriminator loss (LD) and not the power network. As we see in Table 3, the TPR does not suffer from this increased loss for the power network.
For our final hyper-parameter evaluation, we investigate the effect of varying the hinting probability, p, which determines the probability with which Bi = 1. We note that this hyper-parameter is used to trade-off between speed of learning and optimality of the learned solution. A low probability makes for fast convergence, but suboptimal convergence whereas a high probability makes for slow convergence (but a more optimal solution). We chose 0.9 to balance this, following the implementation of [38]. To demonstrate this trade-off we vary this hyper-parameter (from 0 to 0.9). The results below are for the auto-regressive Uniform distribution with Y-Logit and FDR set to 10%.

heading: RESULTS WITH FDR THRESHOLD 5%
text: In all of the synthetic experiments above, the FDR threshold is set to 10% (which is the level most thoroughly investigated by [3]). In this final experiment, we set the FDR threshold to 5% (which aligns with our real-world experiment) and verify that KnockoffGAN is capable of controlling the FDR at this level and still outperforms the original knockoff framework. We conduct the experiments on (1) an auto-regressive distribution with Gaussian marginal distributions, (2) an auto-regressive distribution with Uniform marginal distributions, and (3) a 4-mixture Gaussian mixture model.

heading: STATISTICS OF REAL-WORLD BIOBANK DATASET
text: To preserve anonymity of the authors, the full details of this dataset will be given upon acceptance of the paper. In Table 1, we provide the basic statistics of the real-world biobank dataset.

Review summaries:
The paper presents a deep-learning-based version of the knockoff method by Candes et al. for FDR control in feature selection problems to avoid assumptions posed on the distribution of features by the original method. In a supervised feature selection setting, the goal of the knockoff framework is to select a set of input features that are statistically associated to an output variable Y, while controling the FDR. The basic idea behind knockoff is to generate artificial input feature vectors, (i.e. knockoffs) that are independent of Y, when conditioned on the real feature vector X, but after swapping arbitrary elements with X, are distributed as X. Sets of associated features and FDR estimates are obtained by contrasting suited feature selction criteria that measure associations of knockoffs and real features with the target Y.

This paper introduces a novel feature selection method by utilizing GAN to learn the distributions. The novelty of this paper is to incorporate two recent works, i.e. knockoff for feature selection and W-GAN for generative models. Compared to the latest knockoff work which requires a known multivariate Gaussian distribution for the feature distribution, the proposed work is able to generate knockoffs for any distribution and without any prior knowledge of it.