Chung and Romano (2013) Exact and Asymptotically Robust Permutation Tests, AOS.

Introduction

  •  X_1, \ldots, X_m \sim P IID.
  •  Y_1, \ldots, Y_n \sim Q IID.  X Y は independent.
  • Let  N = m + n and write

\begin{align*}Z = (Z_1, \ldots, Z_N) = (X_1, \ldots, X_m, Y_1, \ldots, Y_n)\end{align*}

  •  \bar \Omega = \{(P,Q) : P = Q\}.  (P,Q) \in \bar \Omega であるとき, (Z_1, \ldots, Z_N) の同時分布は任意の  (Z_{\pi(1)}, \ldots, Z_{\pi(N)}) の分布に等しい。ただし, \pi(1), \ldots, \pi(N) 1, \ldots, N の permutation. 
  •  \mathbf{G}_N をすべての permutation の集合とする。( |\mathbf{G}_N| = N!
  • 検定統計量  T_{m,n} = T_{m,n}(Z_1, \ldots, Z_N). 原則的にはパワーが最大になるように  T を選びたい。
  •   T_{m,n} = T_{m,n}(Z_{\pi(1)}, \ldots, Z_{\pi(N)}) をすべての permutation について計算し,

\begin{align*}T_{m,n}^{(1)} \le T_{m,n}^{(2)} \le \cdots \le T_{m,n}^{(N!)} \end{align*}

と並べる。

  •  k =N! - [\alpha N !] for a nominal level  0 \lt \alpha \lt 1
  •  \phi(Z) = 1 if  T_{m,n}(Z) \gt T_{m,n}^{(k)}(Z),  \phi(Z) = 0 if  T_{m,n}(Z) \lt T_{m,n}^{(k)}(Z) とおく。Tieの場合は論文参照。
  • このとき,for any  (P,Q) \in \bar \Omega,

\begin{align*}E_{P,Q}[\phi(Z)] = \alpha \end{align*}

  • また, \hat R_{m,n}^T(t) = \frac{1}{N!}\sum_{\pi \in \mathbf{G}_N}I(T_{m,n}(Z_{\pi(1)}, \ldots, Z_{\pi(N)}) \le t) とおく。

Permutation test: Reject  H_0 if  T_{m,n}(Z) \gt T_{m,n}^{(k)}(Z) or  T_{m,n}(Z) \gt 1 - \alpha \text{ quantile of } \hat R_{m,n}^T(\cdot)

  • 帰無仮説 H_0: (P,Q) \in \Omega_0 where  \Omega_0 \subset \bar \Omegaについて,permutation testは必ずちょうど \alphaの棄却率を有する = Exactness
  •  \Omega_0 \bar \Omega より大きい時に問題が生じる。なぜなら帰無仮説のもとでpermutateしたデータの分布が元のデータの分布と異なるから。
  • たとえば \Omega_0 = \{(P,Q) : \mu(P) = \mu(Q)\}。検定統計量としては  T_{m,n} = \sqrt{N}(\bar X_m - \bar Y_n)。このとき,permutation testは \Omega_0 = \{(P,Q) : \mu(P) = \mu(Q), P \neq Q\}については検定力がない。
  • Neuhaus (1993): by proper studentization of a test statistic, the permutation test can result in asymptotically valid inference even when the underlying distributions are not the same.
  • The goal of this paper: we would like to retain the exactness property when  P = Q, and also have the asymptotic rejection probability be  \alpha for the more general null hypothesis specifying the parameter.

Robust studentized two-sample test

  •  \theta(\cdot) を real-valued parameterとする。興味のある帰無仮説 H_0: \theta(P) = \theta(Q)
  • Assume that  \hat \theta_m = \hat \theta_m(X_1, \ldots, X_m) and  \hat \theta_n = \hat \theta_n(Y_1, \ldots, Y_n) satisfy

\begin{align*}m^{1/2}(\hat \theta_m - \theta(P)) & = \frac{1}{\sqrt{m}}\sum_{i = 1}^m f_P(X_i) + o_P(1) \\ n^{1/2}(\hat \theta_n - \theta(Q)) & = \frac{1}{\sqrt{n}}\sum_{j = 1}^n f_Q(Y_j) + o_P(1) \end{align*}

  • また,これと同様のasymptotic linearityが混合分布  \bar P = p P + q Q からのIIDサンプルについても成り立つと仮定する。

Theorem 2.1: 帰無仮説 H_0: \theta(P) = \theta(Q), 検定統計量を  T_{m,n} = N^{1/2}[\hat \theta_m(X_1, \ldots, X_m) - \hat \theta_n(Y_1, \ldots, Y_n)] とする。このとき, T_{m.n} の permutation distribution について, \sup_t | \hat R_{m,n}^T(t) - \Phi(t/\tau(\bar P))| \overset{p}{\to} 0 が成立する。ただし, \tau^2(\bar P) = \sigma_2(\bar P)/(p(1-p)),  p = \lim n/N.

  • Remark:  T_{m,n} の真の漸近分布は,平均  0 ・ 分散  \frac{1}{p}\sigma^2(P) + \frac{1}{1-p}\sigma^2(Q)正規分布。これは   \hat R_{m,n}^T(t) の極限と異なる。

Theorem 2.2:  \sigma(P) \sigma(Q) の一致推定量として  \hat \sigma_m \hat \sigma_n が得られるとする。また, V_{m,n} = \sqrt{\frac{N}{m}\hat \sigma_m^2 + \frac{N}{n}\hat \sigma^2_n } とする。このとき, S_{m,n} = T_{m,n}/V_{m,n} の permutation distribution  R_{m,n}^S(\cdot) について, \sup_t | \hat R_{m,n}^S(t) - \Phi(t)| \overset{p}{\to} 0 が成立する。

  • k-sampleのケース:Theorem 3.1 参照。

Li, Cai & Li (2020) Transfer learning for high-dimensional linear regression: prediction, estimation, and minimax optimality, arXiv.

Introduction

  • Target model:  y_i^{(0)} = x_i^{(0)\top}\beta + \epsilon_i^{(0)},  i = 1, \ldots, n_0.
  •  \beta \in \mathbb{R}^p,  p can be larger than  n_0.  s: number of nonzero elements of  \beta.  s is much smaller than  p.
  • Auxiliary models:  y_i^{(k)} = x_i^{(k)\top}\omega^{(k)} + \epsilon_i^{(k)},  i = 1, \ldots, n_k,  k = 1, \ldots, K.
  •  \beta \omega^{(k)} は一般的に異なる.しかしもし両者が近い値をとるならば,target model をより効率的に推定できるかもしれない.
  • Let  \delta^{(k)} := \beta - \omega^{(k)}. "Informative" auxiliary samples を次の集合で定義する:

\begin{align*} \mathcal{A}_q := \{k : ||\delta^{(k)}||_q \le h \}, \;\; \text{for} \;\; q \in [0,1] \end{align*}

  •  h が小さいほど  \mathcal{A}_q は informative.  \mathcal{A}_q空集合であることも許容

Estimation with known informative auxiliary samples

  • Let  \mathcal{A} \equiv \mathcal{A}_1.
Oracle trans-Lasso algorithm
  1.  \hat{\omega}^\mathcal{A} := \arg\min_\omega\frac{1}{2(n_\mathcal{A} + n_0)}\sum_{k \in \mathcal{A} \cup\{0\}} ||y^{(k)} - X^{(k)}\omega||^2_2 + \lambda_\omega ||\omega||_1 for some  \lambda_\omega = c\sqrt{\log p/(n_\mathcal{A} + n_0)}
  2. Let  \hat \beta := \hat{\omega}^\mathcal{A} + \hat{\delta}^\mathcal{A}, where  \hat{\delta}^\mathcal{A} := \arg\min_\delta\frac{1}{2 n_0} ||y^{(0)} - X^{(0)}(\hat{\omega}^\mathcal{A} + \delta)||^2_2 + \lambda_\delta ||\delta||_1 for some  \lambda_\delta = c\sqrt{\log p/(n_\mathcal{A} + n_0)}.
  •  \hat{\omega}^\mathcal{A} の probability limit を  \omega^\mathcal{A} = \beta + \delta^\mathcal{A} とする.また, ||\delta^\mathcal{A}||_1 \le h である.
  • Parameter space:  \Theta_q(s,h) := \{(\beta, \delta^{(1)}, \ldots , \delta^{(K)}) : ||\beta||_0 \le s, \; \max_{k \in \mathcal{A}_q}||\delta^{(k)}||_q \le h\}
  •  k \in \mathcal{A}\cup \{0\} について, X_i^{(k)} \sim i.i.d. Gaussian +  \epsilon_i^{(k)} \sim i.i.d. sub-Gaussian など仮定し,それらの下で以下を得る:

Theorem 1.

 \sup_{\beta \in \Theta_1(s,h)} \frac{1}{n_0}||X^{(0)}(\hat \beta - \beta)||^2_2 \lor ||\hat \beta - \beta||_2^2 = O_P\left( \frac{s \log p}{n_\mathcal{A} + n_0} + \frac{s \log p}{n_0} \land h \sqrt{\frac{\log p}{n_0}} \land h^2\right)

  •  \mathcal{A}空集合の時は右辺は  O_P(s \log p/n_0).これは通常のLassoと同じレート.
  •  h \lt\lt s\sqrt{\log p/n_0} ならば転移学習によって効率性が改善する.つまり, \delta のほうが  \beta よりはるかに sparse である必要.

Theorem 2.

 \inf_{\hat \beta}\sup_{\beta \in \Theta_1(s,h)}\Pr\left( ||\hat \beta - \beta||_2^2 \ge c_1 \frac{s \log p}{n_\mathcal{A} + n_0} + c_2 \frac{s \log p}{n_0} \land h \sqrt{\frac{\log p}{n_0}} \land h^2\right) \ge \frac{1}{2}

  • Theorem 1 + Theorem 2 => Oracle trans-Lasso は minimax rate optimal.

Unknown Set of Informative Auxiliary Samples

  • Let  \hat Q(\mathcal{I}, b) := \sum_{i \in \mathcal{I}}||y_i^{(0)} - x_i^{(0)\top}b||_2^2 denote the sum of squared prediction error under  b, and  \Lambda^{L + 1} := \{\nu \in \mathbb{R}^{L + 1}: \nu_l \ge 0, \: \sum_{l=0}^L \nu_l = 1\} denote an  L-dimensional simplex.
Oracle trans-Lasso algorithm
  1. Let  \mathcal{I} be a random subset of  \{1, \ldots, n\}.
  2. Construct  L + 1 candidate sets of  \mathcal{A},  \{\hat G_0, \hat G_1, \ldots , \hat G_L\} with  \hat G_0 = \emptyset.
  3.  \hat G_l \{X_\mathcal{I}^{(0)}, y_\mathcal{I}^{(0)}\} を用いて trans-Lasso を回す.その結果を  \hat \beta (\hat G_l) とする.
  4. Compute  \hat \theta := \arg\min_{\theta \in \Lambda^{L + 1}} \hat Q(\mathcal{I}^c, \sum_{l=0}^L \theta_l \hat \beta(\hat G_l)) + \sum_{l = 0}^L \theta_l \hat Q(\mathcal{I}^c, \hat \beta(\hat G_l)) + \frac{2 \lambda_\theta}{n_0} \sum_{l=0}^L \theta_l \log (\theta_l). Output  \hat \beta^{\hat \theta} := \sum_{l = 0}^L \hat \theta_l \hat \beta(\hat G_l).

 

  •  \hat G_0空集合にしておくことで,auxiliary sample が完全に uninformative でも結果が(そこまで)悪化しないように保険.
  •  \Pr(\hat G_l \subseteq \mathcal{A}, \;\; \text{for some}\;\; 1\le l \le L) \to 1 ならば Oracle trans-Lasso に近い結果が期待できる.
  •  \hat G_l として  2^K すべての組み合わせを考えるのは計算が大変(理論的にはOK??).
  • どうやって  \hat G_l を決めればいい?
    •  \mathcal{A} の構成により, \{\delta^{(k)}\}_{k \in \mathcal{A}} \{\delta^{(k)}\}_{k \in \mathcal{A}^c} では前者のほうがより sparse.
    •  \delta^{(k)} が sparse な  k-th sample から  \hat G_l に追加すればよい.ただし  \delta^{(k)} は未知なので  \{X^{(k)}, y^{(k)}\} \{X_\mathcal{I}^{(0)}, y_\mathcal{I}^{(0)}\} を用いて推定.具体的な手順は論文参照.

諸々の仮定の下で以下を得る:

Theorem 3.

 \frac{1}{|\mathcal{I}^c|}||X_{\mathcal{I}^c}^{(0)}(\hat \beta^{\hat \theta} - \beta)||^2_2 \lor ||\hat \beta^{\hat \theta} - \beta||_2^2 = O_P\left( \frac{s \log p}{n_{\mathcal{A}^o} + n_0} + \frac{s \log p}{n_0} \land h \sqrt{\frac{\log p}{n_0}} \land h^2 + \frac{\log K}{n_0}\right)

Buhlmann and van de Geer (2011) Section 6.2 Least squares and the Lasso 2/2

6.2.2 続き

  • Let  S \subset \{1, \ldots, p \},  \beta_{j,S} := \beta_j 1\{j \in S\},  \beta_{j,S^c} := \beta_j 1\{j \notin S\}.
  • 明らかに  \beta = \beta_S + \beta_{S^c}.

Lemma 6.3: On  \mathcal{T} with  \lambda \ge 2\lambda_0,  2|| \mathbf{X}(\hat \beta - \beta^0) ||_2^2/n + \lambda || \hat \beta_{S_0^c} ||_1 \le 3 \lambda || \hat \beta_{S_0} -\beta_{S_0}^0 ||_1.

Proof.

Lemma 6.1 と  \lambda \ge 2\lambda_0 から,

\begin{align*}2||\mathbf{X}(\hat \beta - \beta^0) ||_2^2/n + 2\lambda || \hat \beta ||_1 \le \lambda ||\hat \beta - \beta^0|| + 2 \lambda || \beta^0 ||_1\end{align*}

また,三角不等式から

\begin{align*} || \hat \beta ||_1 = || \hat \beta_{S_0} ||_1 + || \hat \beta_{S_0^c} ||_1 \ge || \beta_{S_0}^0 ||_1 - || \hat \beta_{S_0} - \beta_{S_0}^0 ||_1 + || \hat \beta_{S_0^c} ||_1 \end{align*}

さらに, || \beta^0 ||_1 = || \beta^0_{S_0} ||_1,  ||\hat \beta - \beta^0||_1 = ||\hat \beta_{S_0} - \beta^0_{S_0}||_1 + ||\hat \beta_{S_0^c}||_1 から結果を得る.

  • Cauchy-Schwarz inequality から, ||\hat \beta_{S_0} - \beta_{S_0}^0||_1 \le \sqrt{s_0} ||\hat \beta_{S_0} - \beta_{S_0}^0||_2
  • Recall that  \hat \Sigma := \mathbf{X}^\top \mathbf{X}/n,  ||\mathbf{X}(\hat \beta - \beta^0) ||_2^2/n = (\hat \beta - \beta^0)^\top \hat \Sigma (\hat \beta - \beta^0).
  • ここで,ある  \phi_0 \gt 0 について  ||\hat \beta_{S_0} - \beta_{S_0}^0||_2^2 \le (\hat \beta - \beta^0)^\top \hat \Sigma (\hat \beta - \beta^0)/\phi_0^2 ならば,不等式をさらに進めることができる.
  • 一方, \hat \beta はランダムなので上の不等式を満たす  \phi_0 を直接仮定することはできない. \beta のある集合上で一様に不等式を満たす  \phi_0 を考える.
  • Lemma 6.3 から, \mathcal{T} の上で  ||\hat \beta_{S_0^c}||_1 \le 3 ||\hat \beta_{S_0} - \beta_{S_0}^0||_1 を満たす.

Compatibility condition: ある  \phi_0 \gt 0 ||\beta_{S_0^c}||_1 \le 3 ||\beta_{S_0}||_1 について, ||\beta_{S_0}||_1^2 \le \left( \beta^\top \hat \Sigma \beta \right) s_0 / \phi_0^2.

  • Compatibility condition はどのようなときに満たされる?
    •  ||\beta_{S_0}||_1^2 \le s_0 ||\beta_{S_0}||_2^2 から, \hat \Sigma の最小固有値が十分に大きければOK.
    •  ||\beta_{S_0^c}||_1 \le 3 ||\beta_{S_0}||_1 の制約があるので,直接固有値に仮定を入れるよりは弱い.
    •  S_0 は観測できないので,仮定を直接検証することはできない.
    • Bickel et al. (2009) では restricted eigenvalue assumption と呼んでいる.

Theorem 6.1: Compatibility condition +  \mathcal{T} +  \lambda \ge 2\lambda_0  \Longrightarrow.

\begin{align*}||\mathbf{X}(\hat \beta - \beta^0) ||_2^2/n \le 4 \lambda^2 s_0 / \phi_0^2, \quad ||\hat \beta - \beta^0||_1 \le 4 \lambda s_0 / \phi_0^2 \end{align*}

Proof. 

Lemma 6.3 から

\begin{align*} 2||\mathbf{X}(\hat \beta - \beta^0) ||_2^2/n + \lambda||\hat \beta - \beta^0||_1 & = 2||\mathbf{X}(\hat \beta - \beta^0) ||_2^2/n + \lambda||\hat \beta_{S_0} - \beta^0_{S_0}||_1 + \lambda ||\hat \beta_{S_0^c}||_1 \\ & = 4 \lambda ||\hat \beta_{S_0} - \beta^0_{S_0}||_1 \\ & \le 4 \lambda \sqrt{s_0} ||\hat \beta_{S_0} - \beta^0_{S_0}||_2 \\ & \le 4 \lambda \sqrt{s_0} ||\mathbf{X}(\hat \beta - \beta^0)||_2/(\sqrt{n}\phi_0) \\ & \le  ||\mathbf{X}(\hat \beta - \beta^0)||_2^2/n + 4 \lambda^2 s_0/\phi_0^2 \end{align*}

最後の不等式では  4uv \le u^2 + 4v^2 を使った.

Buhlmann and van de Geer (2011) Section 6.2 Least squares and the Lasso 1/2

6.2.1 Introduction

  •  Y_i = \sum_{j = 1}^p \beta_j X_i^{(j)} + \epsilon_i,  i = 1, \ldots, n. In matrix notation,  \mathbf{Y} = \mathbf{X}\beta + \epsilon.
  •  \mathbf{X} は固定, \epsilon は i.i.d.  N(0,\sigma^2) とする.

For the moment,  p \le n とおく.

  • Least squares estimator:  \hat b := (\mathbf{X}^\top \mathbf{X})^{-1} \mathbf{X} \mathbf{Y}
  • このとき, || \mathbf{X}(\hat b - \beta^0) ||_2^2/\sigma^2 = \epsilon^\top \mathbf{X} (\mathbf{X}^\top \mathbf{X})^{-1}  \mathbf{X}^\top \epsilon/\sigma^2 \sim \chi_p^2
  • とくに, \frac{E || \mathbf{X}(\hat b - \beta^0) ||_2^2}{n} = \frac{\sigma^2}{n}p もわかる (as  trace[\mathbf{X} (\mathbf{X}^\top \mathbf{X})^{-1}  \mathbf{X}^\top] = p).
  • つまり, \mathbf{X}^\top \mathbf{X}/n = I_p と正規化しておけば,各  \beta^0_j \sigma^2/n (全体としては  (\sigma^2 / n) \times p) の精度で推定可能.

以降では, p \gt n と仮定する.

  • Active set:  S_0 := \{j: \beta_j^0 \neq 0\}, Sparsity index:  s_0 = |S_0|.
  •  S_0 が既知かつ  s_0 \le n ならば,上の議論から,推定二乗誤差は  (\sigma^2 / n) \times s_0 になる.
  • The Lasso:

\begin{align*}\hat \beta := \arg\min_\beta \left\{ \frac{|| \mathbf{Y} - \mathbf{X}\beta||_2^2}{n} + \lambda ||\beta||_1\right\}\end{align*}

  •  \lambda をうまく ( \sigma \sqrt{\log p/n} ) のオーダーで選ぶと...

Oracle Inequality:  \frac{|| \mathbf{X}(\hat \beta - \beta^0) ||_2^2}{n} \le \text{cont.} \frac{\sigma^2 \log p}{n} s_0.

つまり, S_0 がわからない場合でも  \log p だけの損失で済む.

6.2.2 The result assuming the truth is linear

Lemma 6.1 (Basic Inequality)

 || \mathbf{X}(\hat \beta - \beta^0) ||_2^2/n + \lambda || \hat \beta ||_1 \le 2 \epsilon^\top \mathbf{X}(\hat \beta - \beta^0)/n + \lambda ||\beta^0||_1.

Proof.

 \hat \beta の定義から,

\begin{align*}|| \mathbf{Y} - \mathbf{X}\hat \beta||_2^2/n + \lambda ||\hat \beta||_1 \le || \mathbf{Y} - \mathbf{X}\beta^0||_2^2/n + \lambda ||\beta^0||_1 \end{align*}

また, \mathbf{Y} - \mathbf{X}\hat \beta = -\mathbf{X}(\hat \beta - \beta^0) + \epsilon なので,これを代入すれば結果を得る.

  •  2 \epsilon^\top \mathbf{X}(\hat \beta - \beta^0)/nempirical process part などと呼ぶ.
  • Empirical process part は

\begin{align*}2|\epsilon^\top \mathbf{X}(\hat \beta - \beta^0)|/n \le 2 \left( \max_j |\epsilon^\top \mathbf{X}^{(j)}| \right) ||\hat \beta - \beta^0 ||_1/n\end{align*}

で bound できる.

  • 通常, \lambda は empirical process part を dominate するように選ぶ.
  •  \mathcal{T} := \left\{ 2 \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \le \lambda_0 \right\} とするとき, \lambda \ge 2\lambda_0 とおくことで  || \mathbf{X}(\hat \beta - \beta^0) ||_2^2/n \lesssim \lambda を目指す.
  • Let  \hat \Sigma := \mathbf{X}^\top \mathbf{X} /n and  \hat{\sigma}^2_j := \hat \Sigma_{j,j} ( j = 1, \ldots, p).

Lemma 6.2

 \hat{\sigma}^2_j = 1 とする.このとき,任意の  t \gt 0 について, \lambda_0 := 2 \sigma \sqrt{\frac{t^2 + 2 \log p}{n}} とおけば,

\begin{align*} \Pr(\mathcal{T}) \ge 1 - 2 \exp[-t^2/2]. \end{align*}

Proof.

仮定から,  V_j := \epsilon^\top \mathbf{X}^{(j)}/\sqrt{n \sigma^2} \sim N(0,1). したがって,

\begin{align*} \Pr(\mathcal{T}) & = \Pr\left(2 \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \gt 2 σ \sqrt{\frac{t^2 + 2 \log p}{n}} \right) \\ & = \Pr\left( \max_j |V_j| \gt \sqrt{t^2 + 2 \log p} \right) \\ & \le 2 p \exp \left[ - \frac{t^2 + 2 \log p}{2}\right] = 2 \exp \left[ - \frac{t^2}{2}\right]. \end{align*}*1

Corollary 6.1 (Consistency of the Lasso)

 \hat{\sigma}^2_j = 1 とする.このとき,任意の  t \gt 0 について, \lambda := 4 \hat \sigma \sqrt{\frac{t^2 + 2 \log p}{n}} とおく.  \hat \sigma \sigma の推定量.このとき,少なくとも  1- \alpha の確率で

\begin{align*} 2||\mathbf{X}(\hat \beta - \beta^0)||^2_2/n \le 3 \lambda ||\beta^0||_1, \end{align*}

ただし  \alpha := 2 \exp[-t^2/2] + \Pr(\hat \sigma \le \sigma)

Proof.

Lemma 6.1 から,

\begin{align*}2 || \mathbf{X}(\hat \beta - \beta^0) ||_2^2/n & \le 4 \epsilon^\top \mathbf{X}(\hat \beta - \beta^0)/n + 2 \lambda  ||\beta^0||_1 - 2 \lambda || \hat \beta ||_1 \end{align*}

ここで,もし  4 \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \le \lambda ならば,三角不等式から

\begin{align*}2 || \mathbf{X}(\hat \beta - \beta^0) ||_2^2/n & \le \lambda || \hat \beta - \beta^0 || + 2 \lambda  ||\beta^0||_1 - 2 \lambda || \hat \beta ||_1 \\ & \le 3\lambda  ||\beta^0||_1\end{align*}

また,

\begin{align*}\Pr(4 \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \le \lambda) & = 1 - \Pr\left(4 \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \gt 4 \hat σ \sqrt{\frac{t^2 + 2 \log p}{n}} \right) \\ & \ge 1 - \Pr\left( \max_j |\epsilon^\top \mathbf{X}^{(j)}|/n \gt \hat σ \sqrt{\frac{t^2 + 2 \log p}{n}},  \hat σ \gt σ \right) - \Pr(\hat σ \le σ) \\ & \ge 1 - 2 \exp[-t^2/2] - \Pr(\hat σ \le σ) \end{align*}

  • この結果から, \lambda \asymp \sqrt{\log p /n} であるとき, ||\beta^0||_1 = o(\sqrt{n/ \log p}) ならば Lasso は一致性を持つ.

*1:Chernoff bound + Standard normality

Gold, Lederer & Tao (2020) Inference for high-dimensional instrumental variables regression, JoE.

Introduction

  • モデル:
    •  \mathbf{y} = \mathbf{X}\beta + \mathbf{u},  \beta = (\beta_j)_{j = 1}^{p_x}
    •  E[\mathbf{u} \mid \mathbf{Z}] = 0
    • Doubly high-dimensional setting:  \mathbf{X},  \mathbf{Z} 両方とも high-dimensional. 
  •  \beta_j について statistical inference したい.
  • van de Geer et al. (2014 AoS)のような二段階のde-biased LASSOを考える.
  • Belloni et al. (2018 arXiv)も似たような設定を考えているが,彼らの方法はNeyman orthogonalityに基づいた別の方法.

Two-stage estimation

  • For each  i = 1, \ldots, n,

\begin{align*}y_i & = x_i^\top \beta + u_i \\x_{ij} & = z_i^\top \alpha^{j} + v_{ij}\end{align*}

  • Exclusion restriction:  E[u_i \mid z_i ] = 0,  E[ v_i \mid z_i ] = 0.
  • In matrix form,

\begin{align*} \mathbf{y} & = \mathbf{X}\beta + \mathbf{u} \\  \mathbf{X} & = \mathbf{Z}\mathbf{A} + \mathbf{V} \equiv   \mathbf{D} + \mathbf{V} \end{align*}

where  \mathbf{X} \in \mathbb{R}^{n \times p_x},  \mathbf{D} = E[\mathbf{X} \mid \mathbf{Z}] \in \mathbb{R}^{n \times p_x},  \mathbf{Z} \in \mathbb{R}^{n \times p_z}, and  \mathbf{A} \in \mathbb{R}^{p_z \times p_x}.

  •  p_x \le p_z for identification.
  • Let  \hat{\Sigma}_z = \mathbf{Z}^\top \mathbf{Z}/n.

Assumption 2.2:  z_i are i.i.d. sub-Gaussian and satisfy  E(z_i) = 0.

Assumption 2.4:  \mathbf{v}^j and  \mathbf{u} are sub-Gaussian.

  • First-stage estimator:  \hat{\alpha}^j \equiv \hat{\alpha}(x^j, \mathbf{Z}),  \hat{\mathbf{A}} = (\hat{\alpha}^1, \ldots, \hat{\alpha}^{p_x}),  \hat{\mathbf{D}} = \mathbf{Z}\hat{\mathbf{A}}.
  • Let  \hat{\Sigma}_d = \hat{\mathbf{D}}^\top \hat{\mathbf{D}}/n.
  • Second-stage estimator:  \hat{\beta} \equiv \hat{\beta}(\mathbf{y}, \hat{\mathbf{D}}).

One-step update

  •  \mathbf{y} = \mathbf{X}\beta + \mathbf{u} をOLS推定する場合

 \tilde{\beta} = \hat{\beta} + \hat{\Theta} \mathbf{X}^\top (\mathbf{y} - \mathbf{X}\hat{\beta})/n, where  \hat{\Theta} is the inverse of  \hat{\Sigma}_x = \mathbf{X}^\top \mathbf{X}/n.

  •  \hat{\Sigma}_x が非特異行列ならば第二項はOLSの性質によりゼロ.一方, p_x \gt n のときは  \hat{\Sigma}_x は可逆でないため,

 \tilde{\beta} = \beta + \hat{\Theta} \mathbf{X}^\top \mathbf{u}/n + \underbrace{(\hat{\Theta}\hat{\Sigma}_x - I_{p_x} )(\beta - \hat{\beta})/n}_{f/\sqrt{n}}.

 \Longrightarrow \sqrt{n}(\tilde{\beta}_j - \beta_j) = \frac{1}{\sqrt{n}} \sum_{i = 1}^n \hat{\theta}_j^\top x_i u_i + f_j.

  •  f o_P(1) であることを確かめればよい.
  •  \hat{\beta} が LASSO のとき, \tilde{\beta} を "desparsified LASSO" や "debiased LASSO" などと呼ぶ.
  • IVモデルに対応する one-step update:  \tilde{\beta} = \hat{\beta} + \hat{\Theta} \hat{\mathbf{D}}^\top (\mathbf{y} - \mathbf{X}\hat{\beta})/n.  \hat{\Theta} \Theta = \Sigma_d^{-1} の推定量
  • Lemma 3.1:  \sqrt{n}(\tilde{\beta} - \beta) = \Theta \mathbf{D}^\top \mathbf{u}/n + \text{remainders}.
  •  \hat{\Theta} をどうやって作る? p_x \gt n のときは  \hat{\Sigma}_d は特異行列.=> Cai et al. (2011 JASA) の CLIME estimator を使う.

Theorem 3.4: いくつかのhigh-level conditionsの下で, \sqrt{n}(\tilde{\beta}_j - \beta_j)/ \omega_j \rightsquigarrow N(0,1).

Two-stage LASSO

  • First-stage estimator:  \hat{\alpha}^j \equiv  \arg\min_{\alpha} || x^j - \mathbf{Z} \alpha||/(2n) + r_j ||\alpha||_1.
  • Second-stage estimator:  \hat{\beta} \equiv  \arg\min_{b} || \mathbf{y} - \hat{\mathbf{D}} b ||/(2n) + r_\beta || b ||_1.
  • 以降は Theorem 3.4 の条件を満たすための tuning parameter  r_j, \; r_\beta の選び方や compatibility condition (see Ch.6, Buhlmann and van de Geer, 2011) などについて.

Lubold, Chandrasekhar & McCormick (2021) Identifying the latent space geometry of network models thorough analysis of curvature, arXiv.

Introduction

  • Hoff, Raftery, and Handcock (2002)のlatent space modelでは,あらかじめ次元と距離の定められた潜在空間を基本に,各個人間のconnectionのdpendenceをモデル化している.
  • この論文では,その潜在空間の幾何構造自体を推定することを考える.
  • 幾何構造の選択がなぜ重要か?
    • 例えば二次元のユークリッド空間を考えると,「お互いの距離が等しい4人グループ」が構成できない.
    • 例えば球面を潜在空間として考えると,二点間の距離の上限があるので,よりネットワークが形成されやすくなる (McCormick and Zheng 2015 JASA)
    • サプライチェーンで離れた上流と下流のネットワークがは形成されにくい.このような場合は双曲空間(hyperbolic space)が合っているはず.
  • グラフ:Graph  G = (V, E),  n = |V|.  G is undirected and unweighted.

モデル:\Pr(G_{ij} = 1 | \nu, z, M^p(\kappa)) = \exp(\nu_i + \nu_j - d_{M^p(\kappa)}(z_i, z_j))

where

リーマン幾何学について

  • Killing-Hopf Theorem: 断面曲率(sectionial curvature)*1が一定=>ユークリッド面,球面,双曲面のいずれか
  •  m \in M^p における接ベクトルの空間を  T_m(M^p) とする.
  • また  g を計量テンソル(リーマン多様体を考えているのでリーマン計量) g_m(u,v) \mapsto \mathbb{R}_{\ge 0} ( u,v \in T_m(M^p)) とする.
  • 断面曲率は以下の式で定義される:

\begin{align}\kappa_m(u,v) = \frac{g_m(R_m(u,v)v, u)}{g_m(u,u)g_m(v,v) - g_m(u,v)^2}\end{align}

このとき, \kappa_m(u,v) は基底  u,v の取り方に依らない.また,断面曲率は一定だと仮定していたので, \kappa_m(u,v) = \kappa としてよい.

  •  M^p(\kappa) = \{x \in \mathbb{R}^{p+1}, Q(x,x) = \kappa^{-1}\} とする.ただし  Q は双線型形式(内積の一般化)
  • このとき, M^p(\kappa) 上の二点  (x,y) の距離は

\begin{align}d_{M^p(\kappa)}(x, y) = \kappa^{-1/2}\cos^{-1}(\kappa Q(x,y))\end{align}

で与えられる.

  • 球面  \mathbb{S}^p =>  Q(x,y) = \sum_{i=1}^{p+1}x_i y_i(標準内積
  • 双曲面  \mathbb{H}^p =>  Q(x,y) = -x_0 y_0 + \sum_{i=1}^{p}x_i y_i
  • ここで, z の分布関数  F_z のサポートが  K 個の異なる点  Z = (z_1, \ldots, z_k) を含むとき,上の式から,距離行列  D が与えられたときに双線型形式を対応させる  K \times K 行列

\begin{align}W_\kappa = \frac{1}{\kappa}\cos(\kappa^{1/2} D) \;\; ( = Q(Z,Z))\end{align}

を定義する.

Proposition 1.1

  •  Z \overset{isom}{\to} \mathbb{R}^p \iff W_0 \text{ is p.s.d.}
  • The smallest p s.t.  Z \overset{isom}{\to} \mathbb{R}^p is  p = \text{rank}(W_0).

Proposition 1.2

  •  Z \overset{isom}{\to} \mathbb{S}^p(\kappa) \iff \text{the signature of } W_\kappa \text{ is } (a, n-a, 0) for some  a\le p+1 and some positive  \kappa.
  •  Z \overset{isom}{\to} \mathbb{H}^p(\kappa) \iff \text{the signature of } W_\kappa \text{ is } (1, n-a-1, a) for some  a\le p and some negative  \kappa.

※The signature of  A:  (a,b,c), where  a= # of positive eigenvalues,  b= # of zero eigenvalues, and  c= # of negative eigenvalues. 

  • 以上の結果から,距離行列を幾何構造と独立に構成できれば,潜在空間をどの空間に埋め込めるかが判明する.=>グラフの構造からlink formationの確率を計算し,それを距離として用いる.

Mapping Graphs to Manifolds via Distance Matrices

  • 純化のため,個人効果が存在しない場合を考える: \nu_i = 0 \; \; \forall i.また, F_z のサポートは  Z = (z_1, \ldots, z_K) であるとする. 
  • 二点  (z_k, z_{k'}) 間の距離を  d_{k,k'} = - \log(p_{k,k'}) で測定する.ただし, p_{k,k'} = \Pr(\text{nodes at } z_k \text{ and } z_{k'} \text{ connect})
  •  z_k に属するノードの集合を  V_k とすれば, p_{k,k'} は以下のように推定できる:

\begin{align} \hat p _{k,k'} = \frac{1}{|V_k| \cdot |V_{k'}|} \sum_{(i,j) \in V_k \times V_{k'}} G_{ij} \end{align}

  • 距離行列は  \hat D = (- \log(\hat p_{k,k'})) として推定.
  • 個人効果が非ゼロで  F_z が連続の一般ケースは省略*2

Estimating the Curvature and Minimal Latent Dimension

  • ユークリッド面の場合は曲率の推定は不要.球面と双曲面の場合についてのみ曲率を推定する.
  • ここで,ある \kappa_0 について  Z\mathbb{S}^p(\kappa_0) に埋め込み可能だとする.すると,Prop 1.2から, W_{\kappa_0} の最小固有値はゼロ: \lambda_1(W_{\kappa_0}) = 0.
  • したがって, \kappa_0 がユニークならば

\begin{align}\kappa_0 = \arg\min_{\kappa > 0} |\lambda_1(\kappa W_\kappa)|\end{align}

である.

  • 実際には  W_\kappa は未知なので, \hat W_\kappa = \frac{1}{\kappa}\cos(\kappa^{1/2} \hat D) で置き換える.
  • 同様に,双曲面の場合はProp 1.2から「上から二番目の固有値がゼロ」であることが条件として得られるので, |\lambda_{K-1}(\kappa W_\kappa)| を最小にすればよい.
  • 曲率推定量の一致性 =>Proposition 3.1
  •  \kappa が推定出来たら  \hat W_{\hat \kappa} のランクから  p を推定する.

以下は幾何構造の選択に関する仮説検定や実証分析.

 

*1:ガウス曲率の多様体版らしい.要確認.

*2:論文p.11参照

Kaji, Manresa & Pouliot (2020) An adversarial approach to structural estimation, arXiv.

Introduction
  • Genaratice adversarial networks (GAN) を使って経済データをシミュレートして構造推定する.
  • GAN = minimax game

 \displaystyle \min_{generator} \max_{discriminator} \text{classification accuracy}

genarator: シミュレーションによって仮想データを生成

discriminator: 仮想データと実データを判別

上の minimax 問題を解くことで,実データに最も近い仮想データを生成する generator を作る.

  • Goodfellow のオリジナルの GAN では generator と discriminator をどちらもニューラルネットワークで作っていたが,ここでは DGP が明示的に仮定されているので

   generator = structural econometric model

Adversarial Estimation Framework
  • 尤度関数を直接計算することは困難だがシミュレーションは可能なモデルを考える.
  • 観測データ: \{X_i\}_{i = 1}^n \overset{iid}{\sim} P_0
  • パラメータ: \theta such that  P_0 = P_{\theta_0}
  • 疑似データ: \{X_i^\theta\}_{i = 1}^m \overset{iid}{\sim} P _\theta.  \{X_i^\theta\} \theta に依存しない潜在変数  \tilde X_i によって決まるとする  X_i^\theta = T_\theta(\tilde X_i).
    • 例:First-price aucion.  X_i = \text{bid},  \tilde X_i = \text{value}
  • Classification:  D: \mathcal{X} \to [0, 1];  \{D(X) = 1\} =  X は真のDGPから生成されたデータと分類, \{D(X) = 0\} =  X は疑似データと分類

The adversarial estimator*1

\begin{align} \hat \theta = \arg\min_{\theta \in \Theta} \max_{D \in \mathcal{D}} \frac{1}{n} \sum_{i = 1}^n \log D(X_i) + \frac{1}{m} \sum_{i = 1}^m \log (1 - D(X_i^\theta))\end{align}

  • 目的関数は  2 \log(1/2) から  0 の間の値に収まる:
    • 下限: P_0 P_\theta が区別できない場合, D(X_i) = D(X_i^\theta) = 1/2
    • 上限: P_0 P_\theta が完全に区別できる場合, D(X_i) =1 かつ  D(X_i^\theta) = 0
  • このminimax 問題の population counterpart は

\begin{align} \min_{\theta \in \Theta} \max_{D \in \mathcal{D}} \mathbb{E}_{P_0}[ \log D(X_i)] + \mathbb{E}_{P_\theta}[ \log (1 - D(X_i^\theta))]\end{align}

と表すことができる.

  • このとき,内側の  D に関する最大化問題は明示的に解くことができて

\begin{align} D_\theta(X) := \frac{p_0(X)}{p_0(X) + p_\theta(X)} \end{align}

となる(Goodfellow et al., 2014, Proposition 1).ただし, p P の density. 

  •  D_\theta(X) D として使うとき,上の population 目的関数は  P_0 P_\thetaJensen-Shannon deivergence に他ならない:

\begin{align} \min_{\theta \in \Theta} \max_{D \in \mathcal{D}} \mathbb{E}_{P_0}[ \log D(X_i)] + \mathbb{E}_{P_\theta}[ \log (1 - D(X_i^\theta))] = \min_{\theta \in \Theta} \text{JS}(P_0 || P_\theta)\end{align}

  • さらに, \theta_0 が識別可能ならば  \theta_0 = \min_{\theta \in \Theta} \text{JS}(P_0 || P_\theta) が成り立つ(Goodfellow et al., 2014, Theorem 1).
  • Example 2 (Logistic discriminator)  \displaystyle D(X) = \frac{\exp(\lambda'X)}{1 + \exp(\lambda' X)}. このとき, \hat \theta について解いてみると, n^{-1}\sum_{i = 1}^n X_i \approx m^{-1}\sum_{i = 1}^m X_i^{\hat \theta} となる.
  • Example 3 (Oracle discriminator)  D_\theta を使って  \theta を推定すると

\begin{align} \hat \theta & = \arg\min_{\theta \in \Theta} \text{sample JS}(P_0 || P_\theta)\\& = \arg\min_{\theta \in \Theta}  \frac{1}{n} \sum_{i = 1}^n \log \frac{p_0(X_i)}{p_0(X_i) + p_\theta(X_i)} + \frac{1}{m} \sum_{i = 1}^m \log \frac{p_\theta(X_i^\theta)}{p_0(X_i^\theta) + p_\theta(X_i^\theta)}\end{align}

このとき, n/m \to 0 の下で  \hat \theta は効率的(MLEに等しい).

Statistical Properties
  • Bauer and Kohler (2019) Definition 2:

定義:Generalized hierarchical interaction model (GHIM)

 d \in \mathbb{N}_0,  d^* \in \{1, \dots , d\},  m(\cdot): \mathbb{R}^d \to \mathbb{R}.

 m = GHIM of order  d^* and level  0

 \iff For all  x \in \mathbb{R}^d,  \exists (a_1, \dots , a_{d^*}) \in \mathbb{R}^d and  f(\cdot): \mathbb{R}^d \to \mathbb{R} such that  m(x) = f(a_1'x, \dots, a_{d^*}'x)

 m = GHIM of order  d^* and level  l + 1

 \iff For all  x \in \mathbb{R}^d,  \exists K \in \mathbb{N}_0,  g_k(\cdot): \mathbb{R}^{d^*} \to \mathbb{R}, and  f_{1k}, \dots , f_{d^* k}: \mathbb{R}^d \to \mathbb{R} such that  f_{1k}, \dots , f_{d^* k} are GHIM of order  l and  m(x) = \sum_{k = 1}^K g_k(f_{1k}(x), \dots, f_{d^*k}(x))

  • Bauer and Kohler (2019) Theorem 1:
 \mathbb{E}[Y \mid X = x] = m(x) = GHIM of order  d^* and finite level  l.  \hat m_n(\cdot) Y X に最小二乗回帰して求めた多層ニューラルネットワークとする*2.このとき, \displaystyle \mathbb{E}|| m(X) - \hat m_n(X)||^2 = O\left( n^{- \frac{2p}{2p + d^*} }\right).  p は smoothness parameter*3.

Assumption 3 + regularity conditions  \Longrightarrow

Theorem 5: Hellinger distance between  p_{\hat \theta} and  p_0 =  O_P(n^{-1/2})

Theorem 6:  \sqrt{n}(\hat \theta - \theta_0) \overset{d}{\to} N(0, I_{\theta_0}^{-1} V I_{\theta_0}^{-1}).

  • 尤度関数が明示的に書けないので,Theorem 6 の分散を直接計算することは難しい.この論文ではブートストラップを用いることを提案.

*1:推定方法の発想は Indirect Inference と同じ.

*2:隠れ層やユニットの数の決め方は論文参照.

*3:いわゆる  p-smooth class.  p = q + s, where  q 回連続微分可能+ q導関数が指数  s のヘルダー連続.