Discussion on "Model Confidence Bounds for Variable Selection" by Yang Li, Yuetian Luo, Davide Ferrari, Xiaonan Hu, and Yichen Qin

TL;DR

This paper comments on the methodology of model confidence bounds for variable selection, discussing its implications and potential improvements in statistical model selection techniques.

Contribution

It provides a critical discussion on the existing model confidence bounds approach, highlighting areas for refinement and future research directions.

Findings

Highlights limitations of current confidence bounds methods

Suggests potential improvements for variable selection accuracy

Encourages further exploration of confidence bounds in model selection

Abstract

This is a comment on the article mentioned in the title.

Equations

$$\hat{r}(\hat{m},\hat{m})$$

$$\underset{\theta \in {\mathbb{R}}^p}{\inf }{P}_{n,\theta ,{\sigma }^2}({\tilde{m}}_L\subseteq m^{\ast }(\theta )\subseteq {\tilde{m}}_U)=o(1)$$

$$P_{n,\theta ,{\sigma }^2}(A)={\sum }_{m\in {\mathcal{M}}_{all}}{\int }_{{\mathbb{R}}^n}{I}_A(y,m)K_n(y,m)d{Q}_{n,\theta ,{\sigma }^2}(y)$$

$$C{P}_{\theta }=(1/200)\sum _{i=1}^{200}I({\hat{m}}_L^{(i)}\subseteq m^{\ast }\subseteq {\hat{m}}_U^{(i)}).$$

Full text

(a) Then $P(\tilde{m}_{L}\subseteq m^{\ast}\subseteq\tilde{m}_{U})=1+o(1)$ and $P(\left|\tilde{m}_{U}\right|-\left|\tilde{m}_{L}\right|=0)=1+o(1)$; in fact, $P(\tilde{m}_{U}=\tilde{m}_{L}=\hat{m}=m^{\ast})=1+o(1)$ holds. Moreover, assuming (A.1) only, this statement trivially continues to hold if $\tilde{m}_{L}$ and $\tilde{m}_{U}$ are both replaced by $\hat{m}$.

(b) Part (a) also holds if we replace $\tilde{m}_{L}$ by $\hat{m}_{L}$ and $\tilde{m}_{U}$ by $\hat{m}_{U}$ (where the number of bootstrap replications $B$ may depend on $n$).

(b) Consider first the case where $B$ does not depend on $n$: Then, under (A.2), the event $A_{n}=\{\hat{m}=\hat{m}^{(1)}=\ldots=\hat{m}^{(B)}\}$, which is the finite intersection $\bigcap_{b=1}^{B}\left\{\hat{m}^{(b)}=\hat{m}\right\}$, has probability converging to one. Note that $\hat{r}(\hat{m},\hat{m})=1$ and $\sum\nolimits_{m\neq\hat{m}}\hat{r}(m,m)=0$ hold on $A_{n}$. Inspection of (2) in Li et al. (2018) hence shows that $\hat{m}_{L}=\hat{m}_{U}=\hat{m}$ holds on $A_{n}$, which, together with (A.1), proves the first claim in Part (b) for fixed $B$. Consider next the case where $B=B_{n}$ diverges to infinity as $n\rightarrow\infty$: Similarly as in Part (a), it suffices to show that $\hat{r}(\hat{m},\hat{m})$ converges to $1$ in probability (noting that $\sum\nolimits_{m\neq\hat{m}}\hat{r}(m,m)\rightarrow 0$ in probability then follows). Note that

[TABLE]

The final term on the right-hand side does not depend on $b$ and coincides with $\tilde{r}(\hat{m},\hat{m})$. It converges to $1$ in probability as shown in the proof of Part (a). Now, conditional on $Y$, the first term on the right-hand side is the average of independent random variables that are centered versions of Bernoulli-distributed variables. Hence, conditional on $Y$, this term has zero expectation and variance bounded by $\left(4B_{n}\right)^{-1}$, which goes to zero as $n\rightarrow\infty$. An application of Chebyshev’s inequality shows that this term hence converges to zero in probability. This completes the proof of the first claim in case $B_{n}$ diverges to infinity. The case of a general sequence $B_{n}$ now follows by a standard subsequence argument. The proof is thus complete as the last claim reduces to the last claim in Part (a).    ∎

Proposition 1 obviously strengthens and extends Theorem 1 of Li et al. (2018). However, this result also raises some obvious concerns: In the limit, the coverage probability of the model confidence set (MCS) defined by $\tilde{m}_{L}$ and $\tilde{m}_{U}$ ($\hat{m}_{L}$ and $\hat{m}_{U}$, respectively) is guaranteed to equal one (and not only to be $\geq 1-\alpha$), and the ‘length’ of this MCS equals zero. In particular, asymptotically this MCS coincides with the set $\left\{\hat{m}\right\}$ and thus provides no more information than the point-estimate $\hat{m}$, which the MCS was intended to improve. Interestingly, the fact that the MCS asymptotically reduces to $\left\{\hat{m}\right\}$ under the assumptions of Proposition 1 seems to have informally been noticed by Li et al. (2018), at least in the context of model selection by the adaptive Lasso; see the discussion in Section 1.5 of the supplementary material of that paper. Apparently, however, this did not raise any red flags. [We also note that – in light of the proof of Part (a) given above – the proof of Theorem 1 in Li et al. (2018) seems to be much too complicated; in fact, it seems to be incorrect: For example, the final inequality in the first display in the proof of Theorem 1 appears to confuse the complement of $\left\{\hat{m}_{L}\subseteq\hat{m}\right\}$ with $\left\{\hat{m}_{L}\supsetneqq\hat{m}\right\}$.]

It turns out that a fixed-parameter asymptotic analysis, as used in Proposition 1 and also in Theorem 1 of Li et al. (2018), is not appropriate here. In fixed-parameter asymptotics, the true parameter $\theta$ is kept fixed while sample size $n$ increases to infinity. At this point, it seems fitting to repeat a warning issued a while ago by Hajek (1971, p. 153):

*“Especially misinformative can be those limit results that are not uniform. Then the limit may exhibit some features that are not even approximately true for any finite $n$ …” *

In fact, in the presence of model selection, fixed-parameter asymptotic results not only can, but actually often will, mislead as has been amply documented in Leeb and Pötscher (2005, 2006a, 2006b, 2008a, 2008b, 2008c); Pötscher (2009); Pötscher and Leeb (2009); Pötscher and Schneider (2009, 2011). This is also the case here: We shall provide a uniform asymptotic analysis which reveals the misleading character of the pointwise asymptotic result above and in Theorem 1 of Li et al. (2018). To this end we have to amend the notation a bit: For given sample size $n$, true regression parameter $\theta\in{\mathbb{R}}^{p}$, and variance parameter $\sigma^{2}\in(0,\infty)$ we shall write $P_{n,\theta,\sigma^{2}}$ for $P$ to emphasize its dependence on the quantities indicated; similarly we shall write $m^{\ast}(\theta)$ to denote the smallest correct model for the parameter $\theta$. [That is, $P_{n,\theta,\sigma^{2}}$ and $m^{\ast}(\theta)$ replace the symbols $P$ and $m^{\ast}$ used in Li et al. (2018).] Conditions (A.1) and (A.2) then become $P_{n,\theta,\sigma^{2}}(\hat{m}\neq m^{\ast}(\theta))=o(1)$ and $P_{n,\theta,\sigma^{2}}(\hat{m}^{(b)}\neq\hat{m})=o(1)$, respectively. Compared to Proposition 1 (and Theorem 1 of Li et al. (2018)) we also add a mild condition on the design matrix $X$.

(a) Then for every $\sigma^{2}\in(0,\infty)$ we have

[TABLE]

as $n\rightarrow\infty$. In other words, the minimal coverage probability of the MCS defined by $\tilde{m}_{L}$ and $\tilde{m}_{U}$ converges to zero. Moreover, assuming (A.1) only, this statement continues to hold if $\tilde{m}_{L}$ and $\tilde{m}_{U}$ are both replaced by $\hat{m}$.

(b) Part (a) also holds if we replace $\tilde{m}_{L}$ by $\hat{m}_{L}$ and $\tilde{m}_{U}$ by $\hat{m}_{U}$ (where the number of bootstrap replications $B$ may depend on $n$).

[TABLE]

where $Q_{n,\theta,\sigma^{2}}$ is the probability measure on $\mathbb{R}^{n}$ defined by the $N(X\theta,\sigma^{2}I_{n})$-distribution, $K_{n}(y,m)$ is the bootstrap distribution on $\mathcal{M}_{all}$ (i.e., corresponds to a Markov-kernel from $\mathbb{R}^{n}$ to $\mathcal{M}_{all}$), and $I_{A}$ is the indicator function of $A$. Now, choose an arbitrary constant $\gamma>0$ and a vector $\theta^{(0)}\in{\mathbb{R}}^{p}$ such that at least one coordinate of $\theta^{(0)}$ equals zero. By Proposition 1, we have $P_{n,\theta^{(0)},\sigma^{2}}(\tilde{m}_{L}=\tilde{m}_{U}=m^{\ast}(\theta^{(0)}))=1+o(1)$. For each $n$, choose $\theta^{(n)}\in{\mathbb{R}}^{p}$ so that $m^{\ast}(\theta^{(n)})\neq m^{\ast}(\theta^{(0)})$ and so that $\|\theta^{(n)}-\theta^{(0)}\|\leq\gamma/\sqrt{n}$. (E.g., obtain $\theta^{(n)}$ from $\theta^{(0)}$ by adding $\gamma/\sqrt{n}$ to one of the coordinates of $\theta^{(0)}$ that are equal to zero.) Because the sequences of measures $P_{n,\theta^{(0)},\sigma^{2}}$ and $P_{n,\theta^{(n)},\sigma^{2}}$ are mutually contiguous (see below), it follows that $P_{n,\theta^{(n)},\sigma^{2}}(\tilde{m}_{L}=\tilde{m}_{U}=m^{\ast}(\theta^{(0)}))=1+o(1)$. By construction, we have $m^{\ast}(\theta^{(0)})\neq m^{\ast}(\theta^{(n)})$, and hence $P_{n,\theta^{(n)},\sigma^{2}}(\tilde{m}_{L}\subseteq m^{\ast}(\theta^{(n)})\subseteq\tilde{m}_{U})=o(1)$. [Mutual contiguity of $P_{n,\theta^{(0)},\sigma^{2}}$ and $P_{n,\theta^{(n)},\sigma^{2}}$ is seen as follows: $Q_{n,\theta^{(0)},\sigma^{2}}$ and $Q_{n,\theta^{(n)},\sigma^{2}}$ are well-known to be mutually contiguous, see, e.g., Lemma A.1 in Leeb and Pötscher (2006b) and note that this lemma, while given for the case $\sigma^{2}=1$, obviously extends to any $\sigma^{2}>0$. The claim then follows from Lemma 3.6 in Leeb and Pötscher (2006b).]

(b) Completely analogous, where now $\mathcal{M}_{all}$ is replaced by the $B$-fold Cartesian product $(\mathcal{M}_{all})^{B}$, $I_{A}(y,m)$ and $K_{n}(y,m)$ are replaced by $I_{A}(y,m_{1},\ldots,m_{B})$ and $K_{n}(y,m_{1},\ldots,m_{B})$, respectively, and where the latter represents the bootstrap distribution of a bootstrap sample of size $B$.    ∎

Proposition 2 shows that the proposed MCS does not perform as desired in that the worst-case coverage probability tends to zero, and not to the nominal coverage probability $1-\alpha$. It also shows that the fixed-parameter asymptotic setting used in Theorem 1 of Li et al. (2018) and in Proposition 1 paints a misleadingly optimistic picture of the actual performance of the MCS procedure. This is in line with similar findings in the context of inference post-model-selection reported earlier in Leeb and Pötscher (2005, 2006a, 2006b, 2008a, 2008b, 2008c); Pötscher (2009); Pötscher and Leeb (2009); Pötscher and Schneider (2009, 2011). In fact, the final claim in Part (a) of Proposition 2 can be easily read-off from results concerning model selection probabilities already given in these references (for a start see, e.g., Section 2.1 of Leeb and Pötscher (2005)).

The proof of Proposition 2 reveals that this result continues to hold if the infimum over $\theta\in{\mathbb{R}}^{p}$ in (1) is further restricted, even in a sample-size dependent (e.g., shrinking) way, as long as the infimum at sample size $n$ is taken over a set containing a point $\vartheta^{(n)}$ so that the sequence $\vartheta^{(n)}$ has the same properties as the sequence $\theta^{(n)}$ that is used in the proof. Also, under assumptions (A.1) and (A.2), a result similar to Proposition 2 continues to hold in sufficiently smooth parametric models (as long as the contiguity argument used in the proof goes through). This, of course, also covers regression models with stochastic regressors.

Proposition 2 is an asymptotic result. To illustrate its import for finite-sample situations we performed a small simulation study extending some of the results presented in Li et al. (2018). Data was generated as described in Section 5 of Li et al. (2018). In particular, we consider scenario (c) of Section 5.1 of that paper with the only difference that we consider not only one but several values for the parameter vector $\theta$. More precisely, we set $B=1000$, $n=300$, $p=15$, $p^{\ast}=6$, $\sigma^{2}=1$, $\rho=0.5$, $\theta_{1}=\dots,\theta_{p^{\ast}-1}=1$, and vary $\theta_{p^{\ast}}$ through the values $\theta_{p^{\ast}}=0.05,0.1,0.2,0.3,0.4,0.5,1,2$. Note that $m^{\ast}=\{1,\dots,p^{\ast}\}$. The nominal confidence level was set to 90%, i.e., $\alpha=0.1$. The SCAD, Lasso, minimum BIC, and minimum AIC were used as model selectors. We used the residual bootstrap for SCAD, BIC, and AIC, and the modified residual bootstrap for Lasso. Descriptions of the bootstrap algorithms are provided in Section 1.3 of the supplementary material of Li et al. (2018).

For each of the eight different values of $\theta$ listed above and for each of the four model selectors, we generated 200 response vectors $Y^{(i)}$, obtained $(\hat{m}_{L}^{(i)},\hat{m}_{U}^{(i)})$, and computed an estimate for the coverage probability via

[TABLE]

Let $CP_{\ast}$ denote the minimum of the $CP_{\theta}$-values when $\theta$ varies through the eight values mentioned above. This results in an estimate of a (loose) upper bound for the worst-case coverage probability for each of the four model selection procedures considered. [Of course, we could have searched more thoroughly over the parameter space for $\theta$ to get a smaller upper bound for the worst-case coverage probability, but this would require a much higher investment in computational resources (including a second-stage simulation correcting downward bias resulting from searching for the minimum). Searching only over eight values of $\theta$, as we do here, incurs only a, for our purpose, negligible bias. The rough upper bound $CP_{\ast}$ we obtain is good enough to make our point.] We summarize the values of $CP_{\ast}$ for the four procedures in the subsequent table:

[TABLE]

The result for $CP_{\ast}$, when the consistent model selection procedure minimum BIC is used, is in good agreement with Proposition 2, as $CP_{\ast}$ is much smaller than the nominal confidence level $0.9$. [In additional simulations for model selection by minimum BIC using other values of $\theta$ we even have found values of $CP_{\theta}$ as small as 0.06.] Also, for SCAD the observed $CP_{\ast}$ is much smaller than $0.9$. [Following Li et al. (2018) we have chosen the tuning parameter via crossvalidation. Li et al. (2018) do not show that this results in a consistent model selection procedure, and probably it does not. However, we have not investigated this issue. If the so-tuned SCAD is not a consistent model selection procedure, then Proposition 2 strictly speaking does not apply. Nevertheless, the simulations show that the worst-case coverage probability can be way below the nominal one.] For the minimum AIC procedure, as well as for Lasso, the observed $CP_{\ast}$ is about right. Since minimum AIC is not a consistent, but rather a conservative, model selection procedure, Proposition 2 does not apply; probably the same is true for Lasso given the tuning used. [Again the tuning parameter for Lasso has been chosen via crossvalidation. Li et al. (2018) do not show that this results in a consistent model selection procedure, and probably it does not. We have not investigated this issue any further.] Of course, this is not to say that the worst-case coverage probability for both of these procedures is in any way guaranteed to be close to the nominal one; it could also be that a more exhaustive search over the parameter space would have turned up a much lower worst-case coverage probability also for these procedures.

The above discussion shows that the proposed MCS has substantial defects when used with consistent model selection procedures. It leaves open the possibility that the MCS procedure suggested by Li et al. (2018) may have merit for some conservative model selection procedures; however, as no theoretical support for this is presented in Li et al. (2018), this still needs to be investigated. Another possible route for research is to consider alternative asymptotic scenarios where, e.g., the number of parameters increases with sample size. Problems of prediction following model selection have been successfully tackled in such a framework, see Leeb (2008, 2009); Steinberger and Leeb (2018).