Sparse Hanson-Wright Inequality for a Bilinear Form of Sub-Gaussian Variables

Seongoh Park School of Mathematics, Statistics and Data Science, Sungshin Women’s University, Seoul, Korea Xinlei Wang Department of Statistical Science, Southern Methodist University, Dallas, TX, USA Johan Lim¹¹1To whom all correspondence should be addressed. Email: johanlim@snu.ac.kr Department of Statistics, Seoul National University, Seoul, Korea

Abstract

In this paper, we derive a new version of Hanson-Wright inequality for a sparse bilinear form of sub-Gaussian variables. Our results are generalization of previous deviation inequalities that consider either sparse quadratic forms or dense bilinear forms. We apply the new concentration inequality to testing the cross-covariance matrix when data are subject to missing. Using our results, we can find a threshold value of correlations that controls the family-wise error rate. Furthermore, we discuss the multiplicative measurement error case for the bilinear form with a boundedness condition.

Keywords: Concentration inequality; covariance matrix; missing data; measurement error; general missing dependency.

1 Introduction

Let $(Z_{1j},Z_{2j})$ , $j=1,\ldots,n$ , be pairs of (possibly dependent) sub-Gaussian random variables with mean zero. Suppose $\gamma_{ij}$ is a random variable that corresponds to $Z_{ij}$ and is independent with $Z_{ij}$ . We allow (general) dependency between $\gamma_{1j}$ and $\gamma_{2j}$ . We write the random variables by vectors $\boldsymbol{Z}_{i}=(Z_{i1},\ldots,Z_{in})^{\top}$ , $\boldsymbol{\gamma}_{i}=(\gamma_{i1},\ldots,\gamma_{in})^{\top}$ , $i=1,2$ . For a given $n\times n$ non-random matrix $\boldsymbol{A}$ , Rudelson and Vershynin, (2013) dealt with a concentration inequality of a quadratic form $\boldsymbol{Z}_{1}^{\top}\boldsymbol{A}\boldsymbol{Z}_{1}$ , which is known as Hanson-Wright inequality. Recently, this result has been extended in two directions. On the one hand, Zhou, (2019) included Bernoulli random variables $\boldsymbol{\gamma}_{1}$ in the quadratic form, i.e., $(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})$ where $*$ is the Hadamard (element-wise) product of two vectors (or matrices). On the other hand, Park et al., (2021) proved that the sub-Gaussian vectors could be distinct $\boldsymbol{Z}_{1}\neq\boldsymbol{Z}_{2}$ and derived the deviation inequality of a bilinear form $\boldsymbol{Z}_{1}^{\top}\boldsymbol{A}\boldsymbol{Z}_{2}$ .

In this paper, considering the two extensions all together, we aim to analyze the concentration of the bilinear form $(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})$ where $\boldsymbol{\gamma}_{1},\boldsymbol{\gamma}_{2}$ are vectors of Bernoulli random variables that may be dependent on each other. Furthermore, we generalize the Bernoulli variables to bounded random variables, which is the first attempt in the literature. Table 1 compares different types of bilinear forms. One may consider converting the bilinear form into a quadratic form by concatenating $\boldsymbol{Z}_{1},\boldsymbol{Z}_{2}$ , i.e.

(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})=(\bar{\boldsymbol{Z}}*\bar{\boldsymbol{\gamma}})^{\top}\begin{bmatrix}\boldsymbol{0}&\boldsymbol{A}/2\\ \boldsymbol{A}/2&\boldsymbol{0}\end{bmatrix}(\bar{\boldsymbol{Z}}*\bar{\boldsymbol{\gamma}})

where $\bar{\boldsymbol{Z}}=(\boldsymbol{Z}_{1}^{\top},\boldsymbol{Z}_{2}^{\top})^{\top},\bar{\boldsymbol{\gamma}}=(\boldsymbol{\gamma}_{1}^{\top},\boldsymbol{\gamma}_{2}^{\top})^{\top}$ . However, we cannot directly apply the previous results (e.g. Theorem 1 in Zhou, (2019)) because entries in $\bar{\boldsymbol{Z}}$ are not independent.

We apply the new version of Hanson-Wright inequality to estimating the high-dimensional cross-covariance matrix when data are subject to either missingness or measurement errors. We treat the problem as multiple testing of individual elements of the matrix and propose an estimator thresholding their sample estimates. To determine the cutoff value, we make use of the extended Hanson-Wright inequality and thus can control the family-wise error rate at a desirable level.

	Identical ( $\boldsymbol{Z}_{1}=\boldsymbol{Z}_{2}$ )	Distinct ( $\boldsymbol{Z}_{1}\neq\boldsymbol{Z}_{2}$ )
Dense	$\boldsymbol{Z}_{1}^{\top}\boldsymbol{A}\boldsymbol{Z}_{1}$ , (Rudelson and Vershynin, 2013)	$\boldsymbol{Z}_{1}^{\top}\boldsymbol{A}\boldsymbol{Z}_{2}$ , (Park et al., 2021)
Sparse	$(\boldsymbol{Z}_{1}\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{1}\boldsymbol{\gamma}_{1})$ , (Zhou, 2019)	$(\boldsymbol{Z}_{1}\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}\boldsymbol{\gamma}_{2})$

Table 1: Comparison of bilinear forms in the literature and their references.

The paper is organized as follows. In Section 2, we present the main result of $(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})$ , Theorem 1, and apply it to the problem of testing the cross-covariance matrix in Section 3. In Section 4, we conclude this paper with a brief discussion.

2 Main result

We first define the sub-Gaussian (or $\psi_{2}$ -) norm of a random variable $X$ in $\mathbb{R}$ by

||X||_{\psi_{2}}=\sup_{p\geq 1}\dfrac{(\mathbb{E}|X|^{p})^{1/p}}{\sqrt{p}},

and $X$ is called sub-Gaussian if its $\psi_{2}$ -norm is bounded. For a matrix $\boldsymbol{A}$ , its operator norm $||\boldsymbol{A}||_{2}$ is defined by the square-root of the largest eigenvalue of $\boldsymbol{A}^{\top}\boldsymbol{A}$ . $||\boldsymbol{A}||_{F}=\sqrt{\sum_{i,j}a_{ij}^{2}}$ . For a vector $\boldsymbol{x}$ , $||\boldsymbol{x}||_{2}$ is a 2-norm and, $\boldsymbol{D}(\boldsymbol{x})$ is a diagonal matrix having its diagonal entries by $\boldsymbol{x}$ .

We now describe the main theorem of this paper.

Theorem 1.

Let $(Z_{1j},Z_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) random variables satisfying $\mathbb{E}Z_{1j}=\mathbb{E}Z_{2j}=0$ , and $||Z_{ij}||_{\psi_{2}}\leq K_{i}$ , $i=1,2$ and $\forall j$ . Also, suppose $(\gamma_{1j},\gamma_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) Bernoulli random variables such that $\mathbb{E}\gamma_{ij}=\pi_{ij}$ and $\mathbb{E}\gamma_{1j}\gamma_{2j}=\pi_{12,j}$ for $j=1,\ldots,n$ , $i=1,2$ . Assume $Z_{ij}$ and $\gamma_{i^{\prime}j^{\prime}}$ are independent for distinct pairs $(i,j)\neq(i^{\prime},j^{\prime})$ . Then, we have that for every $t>0$

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\\ \qquad\qquad\qquad\leq 2\exp\left\{-c\min\left(\dfrac{t^{2}}{K_{1}^{2}K_{2}^{2}||\boldsymbol{A}||_{F,\pi}^{2}},\dfrac{t}{K_{1}K_{2}||\boldsymbol{A}||_{2}}\right)\right\},\end{array}

for some numerical constant $c>0$ . For a matrix $\boldsymbol{A}=(a_{ij},1\leq i,j\leq n)$ , we define $||\boldsymbol{A}||_{F,\pi}=\sqrt{\sum_{j=1}^{n}\pi_{12,j}a_{jj}^{2}+\sum_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2j}}$ .

Note that Theorem 1 is a combination of the results for diagonal and off-diagonal cases given below.

Lemma 1.

Assume $\boldsymbol{Z}_{1},\boldsymbol{Z}_{2},\boldsymbol{\gamma}_{1},\boldsymbol{\gamma}_{2}$ are defined as in Theorem 1. Let $\boldsymbol{A}_{\rm diag}={\rm{diag}}(a_{11},\ldots,a_{nn})$ be a diagonal matrix. We denote $\mathbb{E}\gamma_{1j}\gamma_{2j}=\pi_{12,j}$ . Then, for $t>0$ , we have

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}_{\rm diag}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}_{\rm diag}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\\ \qquad\qquad\qquad\leq 2\exp\left\{-c\min\left(\dfrac{t^{2}}{K_{1}^{2}K_{2}^{2}\sum_{j=1}^{n}\pi_{12,j}a_{jj}^{2}},\dfrac{t}{K_{1}K_{2}\max\limits_{1\leq j\leq n}|a_{jj}|}\right)\right\},\end{array}

for some numerical constant $c>0$ .

Lemma 2.

Assume $\boldsymbol{Z}_{1},\boldsymbol{Z}_{2},\boldsymbol{\gamma}_{1},\boldsymbol{\gamma}_{2}$ are defined as in Theorem 1. Let $\boldsymbol{A}_{\rm off}$ be a $n\times n$ matrix with its diagonals zero. We denote $\mathbb{E}\gamma_{ij}=\pi_{ij}$ for $i=1,2$ and $j=1,\ldots,n$ . Then, for $t>0$ , we have

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}_{\rm off}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\\ \qquad\qquad\qquad\leq 2\exp\left\{-c\min\left(\dfrac{t^{2}}{K_{1}^{2}K_{2}^{2}\sum_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2j}},\dfrac{t}{K_{1}K_{2}||\boldsymbol{A}_{\rm off}||_{2}}\right)\right\},\end{array}

for some numerical constant $c>0$ .

A complete proof of the two lemmas, in spirit of Zhou, (2019), is provided in Appendix. It is noted that our theorem above can yield Theorem 1.1 in Zhou, (2019) and Lemma 5 in Park et al., (2021) under the same setting of each.

To handle more general case where we do not have information about mean, we provide the concentration inequality for non-centered sub-Gaussian variables.

Theorem 2.

Let $(\tilde{Z}_{1j},\tilde{Z}_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) random variables with $\mathbb{E}\tilde{Z}_{1i}=\mu_{1i}$ , $\mathbb{E}\tilde{Z}_{2j}=\mu_{2j}$ , and $||\tilde{Z}_{ij}-\mu_{ij}||_{\psi_{2}}\leq K_{i}$ , $i=1,2$ and $\forall j$ . Also, suppose $(\gamma_{1j},\gamma_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) Bernoulli random variables such that $\mathbb{E}\gamma_{ij}=\pi_{ij}$ and $\mathbb{E}\gamma_{1j}\gamma_{2j}=\pi_{12,j}$ for $j=1,\ldots,n$ , $i=1,2$ . Assume $\tilde{Z}_{ij}$ and $\gamma_{i^{\prime}j^{\prime}}$ are independent for distinct pairs $(i,j)\neq(i^{\prime},j^{\prime})$ . Then, we have that for every $t>0$

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\leq d\exp\left\{-c\min\left(\dfrac{t^{2}}{E_{1}},\dfrac{t}{E_{2}}\right)\right\},\end{array}

for some numerical constants $c,d>0$ . $E_{1}$ and $E_{2}$ are defined by

\begin{array}[]{lll}E_{1}=\max\Bigg{\{}&K_{1}^{2}K_{2}^{2}||\boldsymbol{A}||_{F,\pi}^{2},&V_{1,\pi}^{2}V_{2,\pi}^{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1})\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2})||_{F}^{2},\\ &V_{1,\pi}^{2}K_{2}^{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1}*\boldsymbol{\mu}_{1})^{1/2}\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\pi}_{2})^{1/2}||_{F}^{2},&V_{2,\pi}^{2}K_{1}^{2}||\boldsymbol{D}(\boldsymbol{\pi}_{1})^{1/2}\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2}*\boldsymbol{\mu}_{2})^{1/2}||_{F}^{2},\\ &K_{2}^{2}||\boldsymbol{A}^{\top}(\boldsymbol{\mu}_{1}*\boldsymbol{\pi}_{1})||_{2}^{2},&K_{1}^{2}||\boldsymbol{A}(\boldsymbol{\mu}_{2}*\boldsymbol{\pi}_{2})||_{2}^{2},\\ &V_{1,\pi}^{2}||\{\boldsymbol{A}(\boldsymbol{\mu}_{2}*\boldsymbol{\pi}_{2})\}*\boldsymbol{\mu}_{1}||_{2}^{2},&V_{2,\pi}^{2}||\{\boldsymbol{A}(\boldsymbol{\mu}_{1}*\boldsymbol{\pi}_{1})\}*\boldsymbol{\mu}_{2}||_{2}^{2}\Bigg{\}},\\[10.00002pt] E_{2}=\max\Bigg{\{}&K_{1}K_{2}||\boldsymbol{A}||_{2},&V_{1,\pi}V_{2,\pi}||\boldsymbol{D}(\boldsymbol{\mu}_{1})\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2})||_{2},\\ &V_{1,\pi}K_{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1})\boldsymbol{A}||_{2},&V_{2,\pi}K_{1}||\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2})||_{2}\Bigg{\}}.\end{array}

For a matrix $\boldsymbol{A}=(a_{ij},1\leq i,j\leq n)$ , we define $||\boldsymbol{A}||_{F,\pi}=\sqrt{\sum_{j=1}^{n}\pi_{12,j}a_{jj}^{2}+\sum_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2j}}$ . Also, we define $V_{1,\pi}=\max\limits_{1\leq i\leq n}\pi_{1i}(1-\pi_{1i})$ and $V_{2,\pi}$ similarly.

The bilinear form of $\tilde{Z}_{ij}$ in Theorem 2 can be decomposed into the bilinear forms and the linear combinations of centered random variables. Then, we can apply either Theorem 1 or Hoeffding’s inequality (Theorem 5 in Appendix) to them. The detail of the proof can be found in Appendix B.

3 Application to testing the cross-covariance matrix

A cross-covariance matrix $\boldsymbol{\Sigma}_{XY}$ between two random vectors $\boldsymbol{X}=(X_{1},\ldots,X_{p})^{\top}$ and $\boldsymbol{Y}=(Y_{1},\ldots,Y_{q})^{\top}$ , with its $(k,\ell)$ -th entry being $\sigma_{k\ell}^{XY}=\text{Cov}(X_{k},Y_{\ell})$ , refers to the off-diagonal block matrix of the covariance matrix of $\boldsymbol{Z}=\big{(}{\boldsymbol{X}}^{\top},{\boldsymbol{Y}}^{\top}\big{)}^{\top}$ . It is often considered a less important part in the covariance matrix of $\boldsymbol{Z}$ , $\boldsymbol{\Sigma}_{ZZ}$ , and tends to be overpenalized by shrinkage estimators favoring an invertible estimate. However, it is a crucial statistical summary in some applications. For example, the study of multi-omics data, which aims to explain molecular variations at different molecular levels, receives much attention due to the public availability of biological big data and the covariation between two different data sources is just as important as that within each data source. The main question here is to find pairs (or positions in the matrix) of $X_{k}$ ’s and $Y_{\ell}$ ’s that present a large degree of association, which can be treated by hypothesis testing:

\mathcal{H}_{0,k\ell}:\sigma_{k\ell}^{XY}=0\quad\text{vs.}\quad\mathcal{H}_{1,k\ell}:\text{not }\mathcal{H}_{0,k\ell},

(1)

for $1\leq k\leq p,1\leq\ell\leq q$ .

Testing the cross-covariance matrix has not been much explored in literature. Cai et al., (2019) directly address the problem of estimating the cross-covariance matrix. However, they vaguely assume $q\ll p$ and consider simultaneous testing of hypotheses

\mathcal{H}_{0k}:\sigma_{k1}^{XY}=\ldots=\sigma_{kq}^{XY}=0,\quad\text{vs.}\quad\mathcal{H}_{1k}:\text{not }\mathcal{H}_{0k}

for $k=1,\ldots,p$ . They build Hotelling’s $T^{2}$ statistics for individual hypotheses and decide which rows of $\boldsymbol{\Sigma}_{XY}$ are not 0. Hence, the sparsity pattern in Cai et al., (2019) is not the same as considered in this paper. Moreover, their method cannot address missing data.

Considering (1) is equivalent to $\mathcal{H}_{0,k\ell}:\rho_{k\ell}^{XY}=0$ where $\rho_{k\ell}^{XY}=\text{Cor}(X_{k},Y_{\ell})$ , one can also analyze the sample correlation coefficient $\gamma_{k\ell}$ . Bailey et al., (2019) analyzed a universal thresholding estimator of $\gamma_{k\ell}$ based on its asymptotic normality. Though they are interested in estimation of a large correlation matrix, not a cross-correlation matrix directly, their method can be applied to estimation of the cross-covariance matrix. For example, the proposed estimator would be a $p\times q$ matrix with its component $\gamma_{k\ell}1_{(|\gamma_{k\ell}|>c(n,p,q))}$ , and they aim to find a cutoff value $c(n,p,q)$ to control the family-wise error rate (FWER). However, again, if data are subject to missing, their method is no longer valid.

Here, we address the multiple testing problem for the cross-covariance matrix when some of the data are missing or measured with errors. We apply the modified Hanson-Wright inequalities (Theorem 1, 2, 3, 4) in order to choose a threshold value that controls FWER. More specifically, we derive concentration results for an appropriate cross-covariance estimator $\hat{\sigma}_{k\ell}^{XY}$ in the following form:

{\rm P}\left[\max\limits_{k,\ell}\big{|}\hat{\sigma}_{k\ell}^{XY}-\sigma_{k\ell}^{XY}\big{|}>c(n,p,q)\,\right]\leq\alpha,\quad 0<\alpha<1

under some regularity conditions. We begin with the simplest case where data are completely observed and walk through more complicated cases later.

3.1 Complete data case

We begin with the complete data case.

Assumption 1.

Assume each components of random vectors $\boldsymbol{X}\in\mathbb{R}^{p}$ and $\boldsymbol{Y}\in\mathbb{R}^{q}$ are sub-Gaussian, i.e., it holds for some $K_{X},K_{Y}>0$

\max_{1\leq k\leq p}||X_{k}-\mathbb{E}X_{k}||_{\psi_{2}}\leq K_{X},\quad\max_{1\leq\ell\leq q}||Y_{\ell}-\mathbb{E}Y_{\ell}||_{\psi_{2}}\leq K_{Y}.

(2)

Let us denote the mean vector and cross-covariance matrix of $\boldsymbol{X}$ and $\boldsymbol{Y}$ as follows:

\begin{array}[]{l}\mathbb{E}\boldsymbol{X}=\boldsymbol{\mu}^{X}=(\mu_{k}^{X},1\leq k\leq p)^{\top},\\ \mathbb{E}\boldsymbol{Y}=\boldsymbol{\mu}^{Y}=(\mu_{\ell}^{Y},1\leq\ell\leq q)^{\top},\\ \text{Cov}(\boldsymbol{X},\boldsymbol{Y})=\boldsymbol{\Sigma}_{XY}=(\sigma_{k\ell}^{XY},1\leq k\leq p,1\leq\ell\leq q).\end{array}

(3)

Suppose we observe $n$ independent samples $\{(\boldsymbol{X}_{i},\boldsymbol{Y}_{i})\}_{i=1}^{n}$ of $(\boldsymbol{X},\boldsymbol{Y})$ in Assumption 1. Then, the cross-covariance $\sigma_{k\ell}^{XY}$ can be estimated by

\begin{array}[]{rcl}s_{k\ell}&=&\dfrac{1}{n-1}\sum\limits_{i=1}^{n}(X_{ik}-\bar{X}_{k})(Y_{i\ell}-\bar{Y}_{\ell}).\end{array}

where $\bar{X}_{k}$ , $\bar{Y}_{\ell}$ are the sample means corresponding to $\mu_{k}^{X}$ , $\mu_{\ell}^{Y}$ . We can analyze the concentration of each component of $\boldsymbol{S}_{XY}=(s_{k\ell},1\leq k\leq p,1\leq\ell\leq q)$ using Theorem 1 as its special case where all $\pi$ ’s are 1. We first notice that

s_{k\ell}=\dfrac{1}{n}\sum\limits_{i=1}^{n}(X_{ik}-\mu_{k}^{X})(Y_{i\ell}-\mu_{\ell}^{Y})-\dfrac{1}{n(n-1)}\sum\limits_{i\neq j}(X_{ik}-\mu_{k}^{X})(Y_{j\ell}-\mu_{\ell}^{Y}),

where $\mu_{k}^{X}=\mathbb{E}X_{ik}$ and $\mu_{\ell}^{Y}=\mathbb{E}Y_{j\ell}$ . Hence, we can rewrite the sample cross-covariance estimator in a matrix-form by

s_{k\ell}=\boldsymbol{Z}_{1(k)}^{\top}\boldsymbol{A}\boldsymbol{Z}_{2(\ell)}

where $\boldsymbol{Z}_{1(k)}=(X_{1k}-\mu_{k}^{X},\ldots,X_{nk}-\mu_{k}^{X})^{\top},\boldsymbol{Z}_{2(\ell)}=(Y_{1\ell}-\mu_{\ell}^{Y},\ldots,Y_{n\ell}-\mu_{\ell}^{Y})^{\top}$ , and $\boldsymbol{A}=\{n(n-1)\}^{-1}(n\boldsymbol{\rm I}-\boldsymbol{1}\boldsymbol{1}^{\top})$ . Note that $||\boldsymbol{A}||_{F}=1/\sqrt{n-1}$ and $||\boldsymbol{A}||_{2}=1/(n-1)$ . Then, by Theorem 1, the element-wise deviation inequality for the sample cross-covariance is

{\rm P}\bigg{[}\big{|}s_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>t\bigg{]}\leq 2\exp\left\{-\dfrac{c_{1}(n-1)t^{2}}{K_{X}^{2}K_{Y}^{2}}\right\},\quad\dfrac{t}{K_{X}K_{Y}}<1,

for some numerical constant $c_{1}>0$ . Putting $t=K_{X}K_{Y}\sqrt{\log(2pq/\alpha)/\{c_{1}(n-1)\}}$ and using the union bounds, we can get

{\rm P}\left[\max\limits_{k,\ell}\big{|}s_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>C_{1}K_{X}K_{Y}\sqrt{\dfrac{\log(pq/\alpha)}{n-1}}\,\right]\leq\alpha,

if $n/\log(pq/\alpha)>d_{1}$ for some numerical constants $C_{1},d_{1}>0$ .

3.2 Missing data case

For the case where data are subject to missing, we introduce assumptions for missing indicators.

Assumption 2.

Each component $\delta_{k}^{X}$ of the indicator vector $\boldsymbol{\delta}^{X}=(\delta_{k}^{X},1\leq k\leq p)^{\top}$ corresponding to $\boldsymbol{X}$ is 1 if $X_{k}$ is observed and 0 otherwise. $\boldsymbol{\delta}^{Y}$ is similarly defined. Their moments are given by

\begin{array}[]{l}\mathbb{E}\boldsymbol{\delta}^{X}=\boldsymbol{\pi}^{X}=(\pi_{k}^{X},1\leq k\leq p)^{\top},\\ \mathbb{E}\boldsymbol{\delta}^{Y}=\boldsymbol{\pi}^{Y}=(\pi_{\ell}^{Y},1\leq\ell\leq q)^{\top},\\ \mathbb{E}\boldsymbol{\delta}^{X}(\boldsymbol{\delta}^{Y})^{\top}=\boldsymbol{\pi}^{XY}=(\pi_{k\ell},1\leq k\leq p,1\leq\ell\leq q).\end{array}

Note that the above assumption does not mention independence between components of the indicator vector, which means it allows $\delta_{k}^{X}$ and $\delta_{\ell}^{X}$ , $k\neq\ell$ , to be dependent with each other. This implies multiple components in different positions can be missing together under some dependency structure. Assumption 2 of Park et al., (2021) is equivalent to this, and they called it the general missing dependency assumption. For the missing mechanism, missing completely at random (MCAR) is assumed, that is, $(\boldsymbol{X},\boldsymbol{Y})$ is independent of $(\boldsymbol{\delta}^{X},\boldsymbol{\delta}^{Y})$ . More generally, MCAR can be stated as follows.

Assumption 3 (Assumption 3 of Park et al., (2021)).

An event that an observation is missing is independent of both observed and unobserved random variables.

Suppose $n$ independent samples are generated from the population model under Assumption 1, 2, and 3. Each sample consists of the observational data $(\boldsymbol{X}_{i},\boldsymbol{Y}_{i})$ and their missing indicators $(\boldsymbol{\delta}_{i}^{X},\boldsymbol{\delta}_{i}^{Y})$ . However, due to missingness, we can only observe $\tilde{X}_{ik}=\delta_{ik}^{X}X_{ik}$ , $\tilde{Y}_{j\ell}=\delta_{j\ell}^{Y}Y_{j\ell}$ , for $i,j=1,\ldots,n$ , $k=1,\ldots,p$ , and $\ell=1,\ldots,q$ . We can easily check that

\mathbb{E}\big{[}\sum\limits_{i=1}^{n}\tilde{X}_{ik}\tilde{Y}_{i\ell}\big{]}=n\pi_{k\ell}(\sigma_{k\ell}^{XY}+\mu_{k}^{X}\mu_{\ell}^{Y}),\quad\mathbb{E}\big{[}\sum\limits_{i\neq j}^{n}\tilde{X}_{ik}\tilde{Y}_{j\ell}\big{]}=n(n-1)\pi_{k}^{X}\pi_{\ell}^{Y}\mu_{k}^{X}\mu_{\ell}^{Y}.

From the above observation, we define an estimator of the cross-covariance as follows, which is unbiased for $\sigma_{k\ell}^{XY}$ .

\tilde{s}_{k\ell}=\dfrac{\sum_{i=1}^{n}\tilde{X}_{ik}\tilde{Y}_{i\ell}}{n\pi_{k\ell}}-\dfrac{\sum_{i\neq j}^{n}\tilde{X}_{ik}\tilde{Y}_{j\ell}}{n(n-1)\pi_{k}^{X}\pi_{\ell}^{Y}}=(\tilde{\boldsymbol{Z}}_{1(k)}*\boldsymbol{\delta}_{1(k)})^{\top}\boldsymbol{A}_{k\ell}(\tilde{\boldsymbol{Z}}_{2(\ell)}*\boldsymbol{\delta}_{2(\ell)}).

(4)

The last representation is a bilinear form of $\tilde{\boldsymbol{Z}}_{1(k)},\tilde{\boldsymbol{Z}}_{2(k)},\boldsymbol{\delta}_{1(k)},\boldsymbol{\delta}_{2(k)}$ , and $\boldsymbol{A}_{k\ell}$ defined as below.

\begin{array}[]{l}\tilde{\boldsymbol{Z}}_{1(k)}=({X}_{1k},\ldots,{X}_{nk})^{\top},\tilde{\boldsymbol{Z}}_{2(\ell)}=({Y}_{1\ell},\ldots,{Y}_{n\ell})^{\top},\\ \boldsymbol{\delta}_{1(k)}=(\delta_{1k}^{X},\ldots,\delta_{nk}^{X})^{\top},\boldsymbol{\delta}_{2(\ell)}=(\delta_{1\ell}^{Y},\ldots,\delta_{n\ell}^{Y})^{\top},\\ \boldsymbol{A}_{k\ell}=\left(\dfrac{1}{n\pi_{k\ell}}+\dfrac{1}{n(n-1)\pi_{k}^{X}\pi_{\ell}^{Y}}\right)\boldsymbol{\rm I}-\dfrac{1}{n(n-1)\pi_{k}^{X}\pi_{\ell}^{Y}}\boldsymbol{1}\boldsymbol{1}^{\top}.\end{array}

Thus, we apply Theorem 2 to $\tilde{s}_{k\ell}$ and get, for $t<E_{1,k\ell}/E_{2,k\ell}$

{\rm P}\bigg{[}\big{|}\tilde{s}_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>t\bigg{]}\leq d\exp\left\{-\dfrac{c_{2}t^{2}}{E_{1,k\ell}}\right\},

where $c_{2},d>0$ are some numerical constants and $E_{1,k\ell},E_{2,k\ell}$ are defined below.

\begin{array}[]{l}E_{1,k\ell}=\max\bigg{[}\max\Big{\{}K_{X}^{2}K_{Y}^{2},K_{X}^{2}|\mu_{\ell}^{Y}|^{2},|\mu_{k}^{X}|^{2}K_{Y}^{2},|\mu_{k}^{X}|^{2}|\mu_{\ell}^{Y}|^{2}\Big{\}}\left(\dfrac{1}{n\pi_{k\ell}^{2}}+\dfrac{1}{n(n-1)(\pi_{k}^{X})^{2}(\pi_{\ell}^{Y})^{2}}\right),\\ \qquad\qquad\max\Big{\{}(K_{X})^{2}|\mu_{\ell}^{Y}|^{2},|\mu_{k}^{X}|^{2}(K_{Y})^{2},|\mu_{k}^{X}|^{4},|\mu_{\ell}^{Y}|^{4}\Big{\}}\dfrac{1}{n}\left(\dfrac{1}{\pi_{k\ell}}-\dfrac{1}{\pi_{k}^{X}\pi_{\ell}^{Y}}\right)^{2}\bigg{]},\\ E_{2,k\ell}=\max\Big{\{}K_{X}K_{Y},K_{X}|\mu_{\ell}^{Y}|,|\mu_{k}^{X}|K_{Y},|\mu_{k}^{X}||\mu_{\ell}^{Y}|\Big{\}}\left(\dfrac{1}{n\pi_{k\ell}}+\dfrac{1}{n(n-1)\pi_{k}^{X}\pi_{\ell}^{Y}}\right)\end{array}

Putting $t=\sqrt{E_{1,k\ell}\log(dpq/\alpha)/c_{2}}$ and using the union bound argument to the indices $(k,\ell)$ , we can get for some numerical constants $C_{2},d_{2}>0$

{\rm P}\left[\max\limits_{k,\ell}\big{|}\tilde{s}_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>C_{2}\sqrt{\log(pq/\alpha)\max_{k,\ell}E_{1,k\ell}}\,\right]\leq\alpha,

if $\sqrt{\log(pq/\alpha)}<d_{2}\min_{k,\ell}\sqrt{E_{1,k\ell}}/E_{2,k\ell}$ . A simple calculation leads to

\begin{array}[]{l}\max\limits_{k,\ell}E_{1,k\ell}\leq\dfrac{\{f_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})\}^{2}}{(n-1)(\pi_{\min}^{J}\wedge(\pi_{\min}^{M})^{2})^{2}},\\ \min\limits_{k,\ell}\sqrt{E_{1,k\ell}}/E_{2,k\ell}\geq g_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})\sqrt{n-1}(\pi_{\min}^{J}\wedge(\pi_{\min}^{M})^{2}),\end{array}

where $\pi_{\min}^{J}=\min_{k,\ell}\pi_{k\ell}$ , $\pi_{\min}^{M}=\min(\min_{k}\pi_{k}^{X},\min_{\ell}\pi_{\ell}^{Y})$ , $\mu_{X}=\max_{k}|\mu_{k}^{X}|$ , $\mu_{Y}=\max_{\ell}|\mu_{\ell}^{Y}|$ , $f_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})=\max\{K_{X}K_{Y},\mu_{X}K_{Y},K_{X}\mu_{Y},\mu_{X}\mu_{Y},\mu_{X}^{2},\mu_{Y}^{2}\}$ , and $g_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})=\min\{1,\mu_{X}/K_{Y},\mu_{Y}/K_{X},\mu_{X}\mu_{Y}/(K_{X}K_{Y})\}$ . The superscript “J” and “M” stand for joint and marginal, respectively. Then, we conclude that for some numerical constants $\tilde{C}_{2},d_{2}>0$

{\rm P}\left[\max\limits_{k,\ell}\big{|}\tilde{s}_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>\tilde{C}_{2}f_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})\sqrt{\dfrac{\log(pq/\alpha)}{(n-1)(\pi_{\min}^{J}\wedge(\pi_{\min}^{M})^{2})}}\,\right]\leq\alpha,

if $\dfrac{n-1}{\log(pq/\alpha)}>\dfrac{d_{2}}{g_{2}(K_{X},K_{Y},\mu_{X},\mu_{Y})\pi_{\min}^{J}\wedge(\pi_{\min}^{M})^{2}}$ holds.

3.3 Measurement error case

The missing data problem is a special case of the multiplicative measurement error case if the multiplicative factor only takes either 1 or 0. Under some boundedness condition on the factors, we can extend the current framework to the multiplicative measurement error case, which is given below.

Theorem 3.

Let $(Z_{1j},Z_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) random variables satisfying $\mathbb{E}Z_{1j}=\mathbb{E}Z_{2j}=0$ , and $||Z_{ij}||_{\psi_{2}}\leq K_{i}$ , $i=1,2$ and $\forall j$ . Also, suppose $(\gamma_{1j},\gamma_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) non-negative random variables such that $\gamma_{ij}\leq B_{ij}$ almost surely for some $B_{ij}>0$ , $i=1,2$ and $\forall j$ . Assume $Z_{ij}$ and $\gamma_{i^{\prime}j^{\prime}}$ are independent for distinct pairs $(i,j)\neq(i^{\prime},j^{\prime})$ . Then, we have that for every $t>0$

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\boldsymbol{Z}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\boldsymbol{Z}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\\ \qquad\qquad\qquad\leq 2\exp\left\{-c\min\left(\dfrac{t^{2}}{K_{1}^{2}K_{2}^{2}||\boldsymbol{D}(\boldsymbol{B}_{1})\boldsymbol{A}\boldsymbol{D}(\boldsymbol{B}_{2})||_{F}^{2}},\dfrac{t}{K_{1}K_{2}||\boldsymbol{A}||_{2}}\right)\right\},\end{array}

for some numerical constant $c>0$ . $\boldsymbol{D}(\boldsymbol{B}_{1})$ is a diagonal matrix with diagonal elements being from $\boldsymbol{B}_{1}=(B_{1i},1\leq i\leq n)^{\top}$ . $\boldsymbol{D}(\boldsymbol{B}_{2})$ is similarly defined.

The proof is straightforward, and thus we outline the idea behind it and omit the detail. In the diagonal part, we need to modify (5) as

\mathbb{E}(\lambda a_{ii}\gamma_{1i}\gamma_{2i}Z_{1i}Z_{2i})^{s}\leq\mathbb{E}\big{(}\lambda|a_{ii}|B_{1i}B_{2i}|Z_{1i}Z_{2i}|\big{)}^{s},\quad\lambda>0.

For the off-diagonal case, a careful investigation into its proof shows that the missing indicators are conditioned in the analysis and the fact they are Bernoulli random variables is not used until Step 4 in Appendix A.3. The result of Lemma 5 (see Step 4 in Appendix A.3) can be extended to the bounded random errors as we can derive

\mathbb{E}\left[\exp\left(\tau\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\leq\exp\left(\tau\sum\limits_{i\neq j}a_{ij}^{2}B_{1i}B_{2i}\right),\quad\tau>0.

Furthermore, we state the result for the case with non-zero means.

Theorem 4.

Let $(\tilde{Z}_{1j},\tilde{Z}_{2j})$ , $j=1,\ldots,n$ , be independent pairs of (possibly correlated) random variables with $\mathbb{E}\tilde{Z}_{1i}=\mu_{1i}$ , $\mathbb{E}\tilde{Z}_{2j}=\mu_{2j}$ , and $||\tilde{Z}_{ij}-\mu_{ij}||_{\psi_{2}}\leq K_{i}$ , $i=1,2$ and $\forall j$ . Also, suppose $(\gamma_{1j},\gamma_{2j})$ be pairs of (possibly correlated) non-negative random variables such that $\gamma_{ij}\leq B_{ij}$ almost surely for some $B_{ij}>0$ , $i=1,2$ and $\forall j$ . Assume $\tilde{Z}_{ij}$ and $\gamma_{i^{\prime}j^{\prime}}$ are independent for distinct pairs $(i,j)\neq(i^{\prime},j^{\prime})$ . Then, we have that for every $t>0$

\begin{array}[]{l}{\rm P}\bigg{[}\big{|}(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})\big{|}>t\bigg{]}\leq d\exp\left\{-c\min\left(\dfrac{t^{2}}{E_{1}},\dfrac{t}{E_{2}}\right)\right\},\end{array}

for some numerical constants $c,d>0$ . $E_{1}$ and $E_{2}$ are defined by

\begin{array}[]{lll}E_{1}=\max\Bigg{\{}&K_{1}^{2}K_{2}^{2}||\boldsymbol{D}(\boldsymbol{B}_{1})\boldsymbol{A}\boldsymbol{D}(\boldsymbol{B}_{2})||_{F}^{2},&\max\limits_{1\leq i\leq n}B_{1i}^{2}K_{2}^{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1}*\boldsymbol{\mu}_{1})^{1/2}\boldsymbol{A}\boldsymbol{D}(\boldsymbol{B}_{2})^{1/2}||_{F}^{2},\\ &\max\limits_{1\leq i\leq n}B_{2i}^{2}K_{1}^{2}||\boldsymbol{D}(\boldsymbol{B}_{1})^{1/2}\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2}*\boldsymbol{\mu}_{2})^{1/2}||_{F}^{2},&K_{2}^{2}||\{\boldsymbol{A}^{\top}(\boldsymbol{\mu}_{1}*\boldsymbol{u}_{1})\}*\boldsymbol{B}_{2}||_{2}^{2},\\ &K_{1}^{2}||\{\boldsymbol{A}(\boldsymbol{\mu}_{2}*\boldsymbol{u}_{2})\}*\boldsymbol{B}_{1}||_{2}^{2},&\max\limits_{1\leq i\leq n}B_{1i}^{2}\max\limits_{1\leq i\leq n}B_{2i}^{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1})\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2})||_{F}^{2}\Bigg{\}},\\[10.00002pt] E_{2}=\max\Bigg{\{}&K_{1}K_{2}||\boldsymbol{A}||_{2},&\max\limits_{1\leq i\leq n}B_{1i}K_{2}||\boldsymbol{D}(\boldsymbol{\mu}_{1})\boldsymbol{A}||_{2},\\ &&\quad\max\limits_{1\leq i\leq n}B_{2i}K_{1}||\boldsymbol{A}\boldsymbol{D}(\boldsymbol{\mu}_{2})||_{2}\Bigg{\}}.\end{array}

The rest of arguments are similar to Section 3.2. Assume we observe $\tilde{\boldsymbol{X}}_{i}=\boldsymbol{\delta}_{i}^{X}*\boldsymbol{X}_{i}$ and $\tilde{\boldsymbol{Y}}_{i}=\boldsymbol{\delta}_{i}^{Y}*\boldsymbol{Y}_{i}$ , $i=1,\ldots,n$ . While $(\boldsymbol{X}_{i},\boldsymbol{Y}_{i})$ is an independent copy of $(\boldsymbol{X},\boldsymbol{Y})$ in Assumption 1, $(\boldsymbol{\delta}_{i}^{X},\boldsymbol{\delta}_{i}^{Y})$ is no longer a vector of binary variables but an independent copy of $(\boldsymbol{\delta}^{X},\boldsymbol{\delta}^{Y})$ in Assumption 4.

Assumption 4.

Each component $\delta_{k}^{X}$ of $\boldsymbol{\delta}^{X}=(\delta_{k}^{X},1\leq k\leq p)^{\top}$ is a measurement error corresponding to $X_{k}$ of $\boldsymbol{X}$ , which is a non-negative random variable satisfying $\delta_{k}^{X}\leq B_{k}^{X}$ almost surely for each $k$ . $\boldsymbol{\delta}^{Y}$ is similarly defined. Their moments are given by

\begin{array}[]{l}\mathbb{E}\boldsymbol{\delta}^{X}=\boldsymbol{u}^{X}=(u_{k}^{X},1\leq k\leq p)^{\top},\\ \mathbb{E}\boldsymbol{\delta}^{Y}=\boldsymbol{u}^{Y}=(u_{\ell}^{Y},1\leq\ell\leq q)^{\top},\\ \mathbb{E}\boldsymbol{\delta}^{X}(\boldsymbol{\delta}^{Y})^{\top}=\boldsymbol{U}^{XY}=(u_{k\ell},1\leq k\leq p,1\leq\ell\leq q).\end{array}

Accordingly, the unbiased estimator for $\sigma_{k\ell}^{XY}$ is

\check{s}_{k\ell}=\dfrac{\sum_{i=1}^{n}\tilde{X}_{ik}\tilde{Y}_{i\ell}}{nu_{k\ell}}-\dfrac{\sum_{i\neq j}^{n}\tilde{X}_{ik}\tilde{Y}_{j\ell}}{n(n-1)u_{k}^{X}u_{\ell}^{Y}}.

In this case, we can derive

{\rm P}\bigg{[}\big{|}\check{s}_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>t\bigg{]}\leq d\exp\left\{-\dfrac{c_{3}t^{2}}{E_{1,k\ell}}\right\},\quad t<E_{1,k\ell}/E_{2,k\ell}

where $c_{3},d>0$ are some numerical constants and $E_{1,k\ell},E_{2,k\ell}$ are defined below.

\begin{array}[]{l}E_{1,k\ell}=\max\bigg{[}\max\Big{\{}K_{X}^{2}K_{Y}^{2}(B_{k}^{X})^{2}(B_{\ell}^{Y})^{2},|\mu_{k}^{X}|^{2}K_{Y}^{2}(B_{k}^{X})^{2}B_{\ell}^{Y},\\ \qquad\qquad\qquad\qquad K_{X}^{2}|\mu_{\ell}^{Y}|^{2}B_{k}^{X}(B_{\ell}^{Y})^{2},|\mu_{k}^{X}|^{2}|\mu_{\ell}^{Y}|^{2}(B_{k}^{X})^{2}(B_{\ell}^{Y})^{2}\Big{\}}\left(\dfrac{1}{nu_{k\ell}^{2}}+\dfrac{1}{n(n-1)(u_{k}^{X})^{2}(u_{\ell}^{Y})^{2}}\right),\\ \qquad\qquad\qquad\qquad\max\Big{\{}K_{X}^{2}|\mu_{\ell}^{Y}|^{2}(B_{k}^{X})^{2}(u_{\ell}^{Y})^{2},|\mu_{k}^{X}|^{2}K_{Y}^{2}(u_{k}^{X})^{2}(B_{\ell}^{Y})^{2}\Big{\}}\dfrac{1}{n}\left(\dfrac{1}{u_{k\ell}}-\dfrac{1}{u_{k}^{X}u_{\ell}^{Y}}\right)^{2}\bigg{]},\\ E_{2,k\ell}=\max\Big{\{}K_{X}K_{Y},K_{X}|\mu_{\ell}^{Y}|B_{k}^{X},|\mu_{k}^{X}|K_{Y}B_{\ell}^{X}\Big{\}}\left(\dfrac{1}{nu_{k\ell}}+\dfrac{1}{n(n-1)u_{k}^{X}u_{\ell}^{Y}}\right).\end{array}

Moreover, we can observe

\begin{array}[]{l}\max\limits_{k,\ell}E_{1,k\ell}\leq\dfrac{\{f_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})\}^{2}}{(n-1)(u_{\min}^{J}\wedge(u_{\min}^{M})^{2})^{2}},\\ \min\limits_{k,\ell}\sqrt{E_{1,k\ell}}/E_{2,k\ell}\geq g_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})\sqrt{n-1}\cdot u_{\min}^{J}\wedge(u_{\min}^{M})^{2},\end{array}

where $u_{\min}^{J}=\min_{k,\ell}u_{k\ell}$ , $u_{\min}^{M}=\min(\min_{k}u_{k}^{X},\min_{\ell}u_{\ell}^{Y})$ , $\mu_{X}=\max_{k}|\mu_{k}^{X}|$ , $\mu_{Y}=\max_{\ell}|\mu_{\ell}^{Y}|$ , $B_{X}=\max_{k}B_{k}^{X}$ , $B_{Y}=\max_{\ell}B_{\ell}^{Y}$ , $g_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})=\min\{1,K_{X}/(B_{X}\mu_{X}),K_{Y}/(B_{Y}\mu_{Y})\}$ , and

\begin{array}[]{l}f_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})=\max\{K_{X}K_{Y}B_{X}B_{Y},\mu_{X}K_{Y}B_{X}B_{Y},K_{X}\mu_{Y}B_{X}B_{Y},\\ \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\mu_{X}\mu_{Y}B_{X}B_{Y},K_{X}\mu_{Y}B_{X}u_{Y},\mu_{X}K_{Y}u_{X}B_{Y},\}.\end{array}

The superscript “J” and “M” stand for joint and marginal, respectively.

Repeating the calculation as in the previous section, we conclude that for some numerical constants $C_{3},d_{3}>0$

{\rm P}\left[\max\limits_{k,\ell}\big{|}\check{s}_{k\ell}-\sigma_{k\ell}^{XY}\big{|}>C_{3}f_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})\sqrt{\dfrac{\log(pq/\alpha)}{(n-1)\cdot(u_{\min}^{J})^{2}\wedge(u_{\min}^{M})^{4}}}\,\right]\leq\alpha,

if $\dfrac{n-1}{\log(pq/\alpha)}\geq\dfrac{d_{3}}{g_{3}(K_{X},K_{Y},\mu_{X},\mu_{Y},B_{X},B_{Y})\cdot u_{\min}^{J}\wedge(u_{\min}^{M})^{2}}$ holds.

4 Discussion

We discuss the generalized Hanson-Wright inequality where the sparse structure and bilinear form are considered for the first time. This extension facilitates an analysis of concentration of the sample (cross-)covariance estimator even when mean parameters are unknown and some of the data are missing. We apply this result to multiple testing of the cross-covariance matrix.

We further consider a measurement error case extended from the missing data case, which is limited to a bounded random variable. It would be interesting to consider more general multiplicative errors as future work; for example, the truncated normal distribution defined on $(0,\infty)$ can be a good example.

References

Bailey et al., (2019) Bailey, N., Pesaran, M., and Smith, L. V. (2019). A multiple testing approach to the regularisation of large sample correlation matrices. Journal of Econometrics, 208(2):507–534.
Cai et al., (2019) Cai, T., Cai, T., Liao, K., and Liu, W. (2019). Large-scale simultaneous testing of cross-covariance matrices with applications to phewas. Statistica Sinica, 29:983–1005.
Park et al., (2021) Park, S., Wang, X., and Lim, J. (2021). Estimating high-dimensional covariance and precision matrices under general missing dependence. Electronic Journal of Statistics, 15(2):4868 – 4915.
Rudelson and Vershynin, (2013) Rudelson, M. and Vershynin, R. (2013). Hanson-wright inequality and sub-gaussian concentration. Electronic Communications in Probability, 18:9 pp.
Vershynin, (2018) Vershynin, R. (2018). High-dimensional probability : an introduction with applications in data science. Cambridge University Press, Cambridge, United Kingdom New York, NY.
Zhou, (2019) Zhou, S. (2019). Sparse hanson–wright inequalities for subgaussian quadratic forms. Bernoulli, 25(3):1603–1639.

Appendix

Appendix A Proof of Theorem 1

A.1 Proof of Lemma 1 (diagonal part)

Define $S_{0}$ by

S_{0}=\sum_{i=1}^{n}a_{ii}\gamma_{1i}\gamma_{2i}X_{1i}Z_{2i}-\mathbb{E}\left[\sum_{i=1}^{n}a_{ii}\gamma_{1i}\gamma_{2i}Z_{1i}Z_{2i}\right].

For $|\lambda|<1/(4eK_{1}K_{2}\max_{i}|a_{ii}|)$ ,

\begin{array}[]{l}\mathbb{E}\exp(\lambda S_{0})\\[8.00003pt] =\dfrac{\mathbb{E}\exp(\lambda\sum_{i=1}^{n}a_{ii}\gamma_{1i}\gamma_{2i}Z_{1i}Z_{2i})}{\exp\left(\lambda\mathbb{E}\left[\sum_{i=1}^{n}a_{ii}\gamma_{1i}\gamma_{2i}Z_{1i}Z_{2i}\right]\right)}\\[8.00003pt] =\prod\limits_{i=1}^{n}\dfrac{\mathbb{E}_{\gamma}\mathbb{E}_{X}\exp(\lambda a_{ii}\gamma_{1i}\gamma_{2i}Z_{1i}Z_{2i})}{\exp\left(\lambda a_{ii}\mathbb{E}\gamma_{1i}\gamma_{2i}\mathbb{E}Z_{1i}Z_{2i}\right)}\\[8.00003pt] =\prod\limits_{i=1}^{n}\dfrac{(1-\pi_{12,i})+\pi_{12,i}\mathbb{E}_{X}\exp(\lambda a_{ii}Z_{1i}Z_{2i})}{\exp\left(\lambda a_{ii}\pi_{12,i}\mathbb{E}Z_{1i}Z_{2i}\right)}\\[8.00003pt] \leq\prod\limits_{i=1}^{n}\dfrac{(1-\pi_{12,i})+\pi_{12,i}(1+\lambda a_{ii}\mathbb{E}Z_{1i}Z_{2i}+16\lambda^{2}a_{ii}^{2}K_{1}^{2}K_{2}^{2})}{\exp\left(\lambda a_{ii}\pi_{12,i}\mathbb{E}Z_{1i}Z_{2i}\right)}\\[8.00003pt] \leq\prod\limits_{i=1}^{n}\dfrac{\exp\left\{\pi_{12,i}(\lambda a_{ii}\mathbb{E}Z_{1i}Z_{2i}+16\lambda^{2}a_{ii}^{2}K_{1}^{2}K_{2}^{2})\right\}}{\exp\left(\lambda a_{ii}\pi_{12,i}\mathbb{E}Z_{1i}Z_{2i}\right)}\\[8.00003pt] \leq\exp\left\{16\lambda^{2}K_{1}^{2}K_{2}^{2}\sum\limits_{i=1}^{n}\pi_{12,i}a_{ii}^{2}\right\}\end{array}

(5)

where the first inequality uses Lemma 3 and the second holds since $1+x\leq e^{x}$ .

Lemma 3.

Assume that $||Z_{1i}||_{\psi_{2}}\leq K_{1},||Z_{2i}||_{\psi_{2}}\leq K_{2}$ for some $K_{1},K_{2}>0$ and $|\lambda|<1/(4eK_{1}K_{2}\max_{i}|a_{ii}|)$ for given $\{a_{ii}\}_{i=1}^{n}\subset\mathbb{R}$ . Then, for any $i$ , we have

\mathbb{E}\exp\left(\lambda a_{ii}Z_{1i}Z_{2i}\right)-1\leq\lambda a_{ii}\mathbb{E}Z_{1i}Z_{2i}+16\lambda^{2}a_{ii}^{2}K_{1}^{2}K_{2}^{2}.

Proof.

First, we define the sub-exponential (or $\psi_{1}$ -) norm of a random variable by

||X||_{\psi_{1}}=\sup_{p\geq 1}\dfrac{(\mathbb{E}|X|^{p})^{1/p}}{p}.

Since the product of sub-Gaussian random variables is sub-exponential, we define $Y_{i}=Z_{1i}Z_{2i}/||Z_{1i}Z_{2i}||_{\psi_{1}}$ . Setting $t_{i}=\lambda a_{ii}||Z_{1i}Z_{2i}||_{\psi_{1}}$ , we get

\begin{array}[]{l}\mathbb{E}\exp\left(t_{i}Y_{i}\right)\\ =1+t_{i}\mathbb{E}Y_{i}+\sum_{s\geq 2}\dfrac{t_{k}^{s}\mathbb{E}Y_{i}^{s}}{s!}\\ \leq 1+t_{i}\mathbb{E}Y_{i}+\sum_{s\geq 2}\dfrac{|t_{k}|^{s}\mathbb{E}|Y_{i}|^{s}}{s!}\\ \leq 1+t_{i}\mathbb{E}Y_{i}+\sum_{s\geq 2}\dfrac{2|t_{k}|^{s}s^{s}}{s!}\\ \leq 1+t_{i}\mathbb{E}Y_{i}+\sum_{s\geq 2}\dfrac{|t_{k}|^{s}e^{s}}{\sqrt{\pi}}\end{array}

In the last two inequalities, we use the following observations. First, note that for the subexponential variable satisfying $||Y_{i}||_{\psi_{1}}=1$ , it holds (see Prop 2.7.1, Vershynin, (2018))

\mathbb{E}|Y_{i}|^{s}\leq 2s^{s},\quad s\geq 1.

Second, we use Stirling’s formula for $s\geq 2$ that

\dfrac{1}{s!}\leq\dfrac{e^{s}}{2s^{s}\sqrt{\pi}}.

If $|\lambda|<1/(4eK_{1}K_{2}\max_{i}|a_{ii}|)$ , then $e|t_{i}|\leq 1/2$ and thus we get

\sum_{s\geq 3}e^{s}|t_{i}|^{s}\leq(e|t_{i}|)^{3}\sum_{s\geq 0}(1/2)^{s}\leq(e|t_{i}|)^{2}.

Using the above, we have

\begin{array}[]{rcl}\mathbb{E}\exp\left(t_{i}Y_{i}\right)&\leq&1+t_{i}\mathbb{E}Y_{i}+2e^{2}|t_{i}|^{2}/\sqrt{\pi}\\ &=&1+\lambda t_{i}\mathbb{E}Z_{1i}Z_{2i}+2e^{2}\lambda^{2}a_{ii}^{2}||Z_{1i}Z_{2i}||_{\psi_{1}}^{2}/\sqrt{\pi}\\ &\leq&1+\lambda t_{i}\mathbb{E}Z_{1i}Z_{2i}+16\lambda^{2}a_{ii}^{2}K_{1}^{2}K_{2}^{2}.\end{array}

∎

Then, for $x>0$ and $0<\lambda<1/(4eK_{1}K_{2}\max_{i}|a_{ii}|)$ , we have

\begin{array}[]{l}{\rm P}(S_{0}>x)\\ ={\rm P}\left(\exp(\lambda S_{0})>\exp(\lambda x)\right)\\ \leq\dfrac{\mathbb{E}\exp(\lambda S_{0})}{\exp(\lambda x)}\\ \leq\exp\left\{-\lambda x+16\lambda^{2}K_{1}^{2}K_{2}^{2}\sum\limits_{i=1}^{n}\pi_{12,i}a_{ii}^{2}\right\}.\end{array}

For the optimal choice of $\lambda$ , that is,

\lambda=\min\left(\dfrac{x}{32eK_{1}^{2}K_{2}^{2}\sum_{i}\pi_{12,i}a_{ii}^{2}},\dfrac{1}{4eK_{1}K_{2}\max_{i}|a_{ii}|}\right)

we can obtain the concentration bound

{\rm P}(S_{0}>x)\leq\exp\left\{-\min\left(\dfrac{x^{2}}{32K_{1}^{2}K_{2}^{2}\sum_{i}\pi_{12,i}a_{ii}^{2}},\dfrac{x}{8eK_{1}K_{2}\max_{i}|a_{ii}|}\right)\right\}.

Considering $-A_{\rm diag}$ instead of $A_{\rm diag}$ in the theorem, we have the same probability bound for ${\rm P}(S_{0}<-x)$ , $x>0$ , which concludes the proof.

A.2 Proof of Lemma 2 (off-diagonal part)

Define $S_{\rm off}$ by

S_{\rm off}=\sum_{1\leq i\neq j\leq n}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j},

whose expectation is zero due to independence across distinct $i$ and $j$ . Then, we claim the lemma below:

Lemma 4.

Assume $||Z_{1i}||_{\psi_{2}}=||Z_{2i}||_{\psi_{2}}=1$ and let $\{a_{ij}\}_{1\leq i\neq j\leq n}\subset\mathbb{R}$ be given. For $|\lambda|<1/(2\sqrt{C_{4}}||A||_{2})$ for some numeric constant $C_{4}>0$ , we have

\mathbb{E}\exp\left(\lambda S_{\rm off}\right)\leq\exp\left(1.44C_{4}\lambda^{2}\sum_{i\neq j}a_{ij}\pi_{1i}\pi_{2j}\right).

The proof is pending until Appendix A.3. Without loss of generality, we can assume $||Z_{1i}||_{\psi_{2}}=||Z_{2i}||_{\psi_{2}}=1$ ; otherwise set $Z_{1i}\leftarrow Z_{1i}/||Z_{1i}||_{\psi_{2}},Z_{2i}\leftarrow Z_{2i}/||Z_{2i}||_{\psi_{2}}$ . Using Lemma 4, we get for $x>0$ and $0<\lambda<1/(2\sqrt{C_{4}}||A||_{2})$

\begin{array}[]{l}{\rm P}(S_{0}>x)\\ ={\rm P}\left(\exp(\lambda S_{0})>\exp(\lambda x)\right)\\ \leq\dfrac{\mathbb{E}\exp(\lambda S_{0})}{\exp(\lambda x)}\\ \leq\exp\left\{-\lambda x+1.44C_{4}\lambda^{2}\sum_{i\neq j}a_{ij}\pi_{1i}\pi_{2j}\right\}.\end{array}

For the optimal choice of $\lambda$ , that is,

\lambda=\min\left(\dfrac{x}{2.88C_{4}\sum_{i\neq j}a_{ij}\pi_{1i}\pi_{2j}},\dfrac{1}{2\sqrt{C_{4}}||A||_{2}}\right)

we can obtain the concentration bound

{\rm P}(S_{0}>x)\leq\exp\left\{-c\min\left(\dfrac{x^{2}}{\sum_{i\neq j}a_{ij}\pi_{1i}\pi_{2j}},\dfrac{x}{||A||_{2}}\right)\right\}.

Considering $-A_{\rm off}$ instead of $A_{\rm off}$ in the theorem, we have the same probability bound for ${\rm P}(S_{0}<-x)$ , $x>0$ , which concludes the proof.

A.3 Proof of Lemma 4

This proof follows the logic of the proof of (9) in Zhou, (2019). We first decouple the off-diagonal sum $S_{\rm off}$ .

Step 1. Decoupling

We introduce Bernoulli variables $\eta=(\eta_{1},\ldots,\eta_{n})^{\top}$ with $\mathbb{E}\eta_{k}=1/2$ for any $k$ . They are independent with each other and also independent of $Z_{1},Z_{2}$ and $\gamma_{1},\gamma_{2}$ . Given $\eta$ , we define $Z_{1}^{\eta}$ by a subvector of $Z_{1}$ at which $\eta_{i}=1$ and $Z_{2}^{1-\eta}$ by a subvector of $Z_{2}$ at which $\eta_{j}=0$ , respectively. Let $\mathbb{E}_{Q}$ be the expectation with respect to a random variable $Q$ where $Q$ can be any of $Z_{1},Z_{2},\gamma_{1},\gamma_{2},\eta$ in our context, or $X,\gamma$ to denote $Z_{1},Z_{2}$ and $\gamma_{1},\gamma_{2}$ together. Define a random variable

S_{\eta}=\sum_{1\leq i,j\leq n}a_{ij}\eta_{i}(1-\eta_{j})\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j}.

Using $S_{\rm off}=4\mathbb{E}_{\gamma}S_{\eta}$ and Jensen’s inequality with $x\mapsto e^{ax}$ for any $a\in\mathbb{R}$ , we get

\mathbb{E}\exp(\lambda S_{\rm off})=\mathbb{E}_{\gamma}\mathbb{E}_{X}\exp(\mathbb{E}_{\eta}4\lambda S_{\eta})\leq\mathbb{E}_{\gamma}\mathbb{E}_{X}\mathbb{E}_{\eta}\exp(4\lambda S_{\eta}).

We condition all the variables except $Z_{2}^{1-\eta}=(Z_{2j},\eta_{j}=0)$ on $\exp(4\lambda S_{\eta})$ and consider its moment generating function denoted by $f(\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta})$

f(\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta})=\mathbb{E}\left(\exp(4\lambda S_{\eta})|\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta}\right).

Note that $S_{\eta}$ can be seen a linear combination of sub-Gaussian variables $Z_{2j}$ , for $j$ such that $\eta_{j}=0$ , that is,

S_{\eta}=\sum_{j:\eta_{j}=0}Z_{2j}\left(\sum_{i:\eta_{i}=1}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}\right),

conditional on all other variables. Thus, the conditional distribution of $S_{\eta}$ is sub-Gaussian satisfying

||S_{\eta}||_{\psi_{2}}\leq C_{0}\sigma_{\eta,\gamma},

where $\sigma_{\eta,\gamma}^{2}=\sum\limits_{j:\eta_{j}=0}\gamma_{2j}\left(\sum\limits_{i:\eta_{i}=1}a_{ij}\gamma_{1i}Z_{1i}\right)^{2}$ . Therefore, we have that for some $C^{\prime}>0$

f(\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta})=\mathbb{E}\left(\exp(4\lambda S_{\eta})|\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta}\right)\leq\exp(C^{\prime}\lambda^{2}\sigma_{\eta,\gamma}^{2}).

(6)

Taking the expectations on both sides, we get

\mathbb{E}_{\gamma}\mathbb{E}_{Z_{1}^{\eta}}f(\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta})\leq\mathbb{E}_{\gamma}\mathbb{E}_{Z_{1}^{\eta}}\exp(C^{\prime}\lambda^{2}\sigma_{\eta,\gamma}^{2})=:\tilde{f}_{\eta}.

Step 2. Reduction to normal random variables

Assume that $\eta$ , $\gamma_{1},\gamma_{2}$ , and $Z_{1}^{\eta}$ are fixed. Let $g=(g_{1},\ldots,g_{n})^{\top}$ be given where $g_{i}$ is i.i.d. from $N(0,1)$ . Replacing $Z_{2j}$ by $g_{j}$ in $S_{\eta}$ , we define a random variable $Z$ by

Z=\sum_{j:\eta_{j}=0}g_{j}\left(\sum_{i:\eta_{i}=1}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}\right).

Due to the property of Gaussian variables, the conditional distribution of $Z$ follows $N(0,\sigma_{\eta,\gamma}^{2})$ . Hence, its conditional moment generating function is given by

\mathbb{E}_{g}\exp(tZ)=\mathbb{E}\left[\exp(tZ)|\eta,\gamma_{1},\gamma_{2},Z_{1}^{\eta}\right]=\exp(t^{2}\sigma_{\eta,\gamma}^{2}/2).

(7)

Now, consider $Z$ as a linear combination of $\{Z_{1i}:\eta_{i}=1\}$ by rewriting it by

Z=\sum_{i:\eta_{i}=1}Z_{1i}\left(\sum_{j:\eta_{j}=0}g_{j}a_{ij}\gamma_{1i}\gamma_{2j}\right),

where $Z_{1i}$ ’s are regarded random variables and others fixed. Then, we can get for some $C_{3}>0$

\mathbb{E}_{Z_{1}^{\eta}}\exp(\sqrt{2C^{\prime}}\lambda Z)=\mathbb{E}\left[\exp(\sqrt{2C^{\prime}}\lambda Z)|g,\eta,\gamma_{1},\gamma_{2}\right]\leq\exp\left(C_{3}\lambda^{2}\sum_{i:\eta_{i}=1}\gamma_{1i}\left(\sum_{j:\eta_{j}=0}g_{j}a_{ij}\gamma_{2j}\right)^{2}\right).

(8)

We now get an upper bound of $\tilde{f}_{\eta}$ from step 1 by using (7) and (8). First, choosing $t=\sqrt{2C^{\prime}}\lambda$ at (7) gives the same term in (6)

\tilde{f}_{\eta}=\mathbb{E}_{\gamma}\mathbb{E}_{Z_{1}^{\eta}}\mathbb{E}_{g}\exp(\sqrt{2C^{\prime}}\lambda Z)=\mathbb{E}_{\gamma}\mathbb{E}_{g}\mathbb{E}_{Z_{1}^{\eta}}\exp(\sqrt{2C^{\prime}}\lambda Z).

Then, we apply (8) to the right-hand side of the above:

\tilde{f}_{\eta}\leq\mathbb{E}_{\gamma}\mathbb{E}_{g}\exp\left(C_{3}\lambda^{2}\sum_{i:\eta_{i}=1}\gamma_{1i}\left(\sum_{j:\eta_{j}=0}g_{j}a_{ij}\gamma_{2j}\right)^{2}\right)=\mathbb{E}\left[\exp\left(C_{3}\lambda^{2}(g^{\top}B_{\eta,\gamma}g)\right)|\eta\right],

(9)

where $B_{\eta,\gamma}$ has its $(k,\ell)$ -th element $\sum\limits_{i:\eta_{i}=1}\gamma_{1i}a_{ik}a_{i\ell}\gamma_{2k}\gamma_{2\ell}$ if $k,\ell\in\{j:\eta_{j}=0\}$ ; $0$ , otherwise. Note that $B_{\eta,\gamma}$ is symmetric.

Step 3. Integrating out the normal random variables

For fixed $\gamma_{1},\gamma_{2},\eta$ , the quadratic form $g^{\top}B_{\eta,\gamma}g$ is identical in distribution with $\sum_{j=1}^{n}s_{j}^{2}g_{j}^{2}$ due to the rotational invariance of $g$ where $s_{j}^{2}$ is the eigenvalue of $B_{\eta,\gamma}$ . Note that the eigenvalues satisfy, for any realization of $\gamma_{1},\gamma_{2},\eta$ ,

\begin{array}[]{rcl}\max\limits_{1\leq j\leq n}s_{j}^{2}&=&||B_{\eta,\gamma}||_{2}\leq||A||_{2}^{2}\\ \sum\limits_{1\leq j\leq n}s_{j}^{2}&=&{\rm{tr}}(B_{\eta,\gamma})=\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}.\end{array}

(10)

Now, using $g_{j}^{2}\sim\chi^{2}_{1}$ and $\mathbb{E}\exp(tg_{j}^{2})\leq(1-2t)^{-1/2}\leq\exp(2t)$ for $t<1/4$ , we get from (9) and (10)

\begin{array}[]{l}\tilde{f}_{\eta}\leq\mathbb{E}\left[\exp\left(C_{3}\lambda^{2}(g^{\top}B_{\eta,\gamma}g)\right)|\eta\right]\\[5.0pt] =\mathbb{E}_{\gamma}\mathbb{E}\left[\exp\left(C_{3}\lambda^{2}(g^{\top}B_{\eta,\gamma}g)\right)|\eta,\gamma_{1},\gamma_{2}\right]\\[5.0pt] =\mathbb{E}_{\gamma}\mathbb{E}\left[\exp\left(\sum\limits_{j=1}^{n}C_{3}\lambda^{2}s_{j}^{2}g_{j}^{2}\right)|\eta,\gamma_{1},\gamma_{2}\right]\\[5.0pt] \leq\mathbb{E}\left[\prod\limits_{j=1}^{n}\exp(2C_{3}\lambda^{2}s_{j}^{2})|\eta\right]\\[5.0pt] =\mathbb{E}\left[\exp\left(2C_{3}\lambda^{2}\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\end{array}

for $\lambda^{2}<1/(4C_{3}||A||_{2}^{2})(<1/(4C_{3}\max_{j}s_{j}^{2}))$ . It is worth noting that the Bernoulli variables $\gamma_{1},\gamma_{2}$ are decoupled by $\eta$ .

Step 4. Integrating out the Bernoulli random variables

We now integrate out $\gamma_{1},\gamma_{2}$ from $\tilde{f}_{\eta}$ by using the following lemma. For $0<\lambda^{2}\leq 1/(8C_{3}||A||_{2}^{2})$ , we have

\tilde{f}_{\eta}\leq\mathbb{E}\left[\exp\left(2C_{3}\lambda^{2}\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\leq\exp\left(2.88C_{3}\lambda^{2}\sum\limits_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2i}\right)

Lemma 5.

For $0<\tau\leq 1/(4||A||_{2}^{2})$ , the conditional expectation satisfies

\mathbb{E}\left[\exp\left(\tau\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\leq\exp\left(1.44\tau\sum\limits_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2i}\right).

Proof.

Due to the independence of $\gamma_{1i},i=1,2,\ldots,n$ , we get

\begin{array}[]{l}\mathbb{E}\left[\exp\left(\tau\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\gamma_{2}^{1-\eta},\eta\right]\\ =\prod\limits_{i:\eta_{i}=1}\mathbb{E}\left[\exp\left(\tau\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\gamma_{2}^{1-\eta},\eta\right]\\ =\prod\limits_{i:\eta_{i}=1}\left[(1-\pi_{1i})+\pi_{1i}\exp\left(\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)\right].\end{array}

Then, we use the local approximation of the exponential function:

e^{x}-1\leq 1.2x,\quad 0\leq x\leq 0.35.

To use it, $\tau$ should satisfies, for given $\eta$ and $\gamma_{2}$ ,

\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}<0.35,

(11)

which leads to

(1-\pi_{1i})+\pi_{1i}\exp\left(\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)\leq 1+1.2\pi_{1i}\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\leq\exp\left(1.2\pi_{1i}\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right),

where we use $1+x\leq e^{x}$ in the last inequality. A sufficient condition for $\tau$ for any realization of $\eta$ and $\gamma_{2}$ in (11) is $\tau\leq 1/(4||A||_{2}^{2})$ since

\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\leq\tau\max_{i}\sum_{j}a_{ij}^{2}\leq\tau||A||_{2}^{2}\leq 0.25.

Hence, for $\tau\leq 1/(4||A||_{2}^{2})$ ,

\begin{array}[]{l}\mathbb{E}\left[\exp\left(\tau\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\\ =\mathbb{E}_{\gamma_{2}^{1-\eta}}\mathbb{E}\left[\exp\left(\tau\sum\limits_{i:\eta_{i}=1}\gamma_{1i}\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\gamma_{2}^{1-\eta},\eta\right]\\ \leq\mathbb{E}\left[\exp\left(\sum\limits_{i:\eta_{i}=1}1.2\pi_{1i}\tau\sum\limits_{j:\eta_{j}=0}a_{ij}^{2}\gamma_{2j}\right)|\eta\right]\\ =\prod\limits_{j:\eta_{j}=0}\mathbb{E}\left[\exp\left(1.2\tau\gamma_{2j}\sum\limits_{i:\eta_{i}=1}\pi_{1i}a_{ij}^{2}\right)|\eta\right]\\ \leq\prod\limits_{j:\eta_{j}=0}\exp\left(1.2^{2}\tau\pi_{2j}\sum\limits_{i:\eta_{i}=1}\pi_{1i}a_{ij}^{2}\right).\end{array}

We apply the similar argument used for $\gamma_{2j},j\in\{j:\eta_{j}=0\}$ to $\gamma_{1i},i\in\{i:\eta_{j}=1\}$ in the last inequality. Note that $\tau\leq 1/(4||A||_{2}^{2})$ is also sufficient since

1.2\tau\sum\limits_{i:\eta_{j}=1}a_{ij}^{2}\pi_{1i}\leq 1.2\tau\max_{j}\sum_{i}a_{ij}^{2}\leq 1.2\tau||A||_{2}^{2}\leq 0.3.

Finally, we observe that

\exp\left(1.2^{2}\tau\sum\limits_{j:\eta_{j}=0}\sum\limits_{i:\eta_{i}=1}a_{ij}^{2}\pi_{1i}\pi_{2j}\right)\leq\exp\left(1.2^{2}\tau\sum\limits_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2j}\right),

which concludes the proof. ∎

Step 5. Combining things together

Assume $|\lambda|<1/(\sqrt{8C_{3}}||A||_{2})$ . Combining all the steps above, we get

\begin{array}[]{l}\mathbb{E}\exp(\lambda S_{\rm off})\\ \leq\mathbb{E}_{\eta}\mathbb{E}_{\gamma}\mathbb{E}_{X}\exp(4\lambda S_{\eta})\\ \leq\mathbb{E}_{\eta}\mathbb{E}_{\gamma}\mathbb{E}_{Z_{1}^{\eta}}\mathbb{E}\left(\exp(4\lambda S_{\eta})|\gamma_{1},\gamma_{2},\eta,Z_{1}^{\eta}\right)\\ \leq\mathbb{E}_{\eta}\tilde{f}_{\eta}\\ \leq\exp\left(2.88C_{3}\lambda^{2}\sum\limits_{i\neq j}a_{ij}^{2}\pi_{1i}\pi_{2i}\right),\end{array}

which proves Lemma 4.

Appendix B Proof of Theorem 2

Let us define the centered sub-Gaussian variable $Z_{ij}=\tilde{Z}_{ij}-\mu_{ij}$ . Then, it is easy to check that the bilinear form is decomposed into 8 terms.

\begin{array}[]{rcl}S&\equiv&(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})-\mathbb{E}(\tilde{\boldsymbol{Z}}_{1}*\boldsymbol{\gamma}_{1})^{\top}\boldsymbol{A}(\tilde{\boldsymbol{Z}}_{2}*\boldsymbol{\gamma}_{2})\\ &=&\sum\limits_{i,j}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j}-\mathbb{E}\left[\sum\limits_{i,j}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j}\right]\\ &&+\sum\limits_{i,j}a_{ij}\mu_{1i}\gamma_{1i}\gamma_{2j}Z_{2j}\\ &&+\sum\limits_{i,j}a_{ij}\mu_{2j}\gamma_{1i}\gamma_{2j}Z_{1i}\\ &&+\sum\limits_{i,j}a_{ij}\mu_{1i}\mu_{2j}\gamma_{1i}\gamma_{2j}-\mathbb{E}\left[\sum\limits_{i,j}a_{ij}\mu_{1i}\mu_{2j}\gamma_{1i}\gamma_{2j}\right]\\ &=&\sum\limits_{i,j}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j}-\mathbb{E}\left[\sum\limits_{i,j}a_{ij}\gamma_{1i}\gamma_{2j}Z_{1i}Z_{2j}\right]\\ &&+\sum\limits_{i,j}a_{ij}\mu_{1i}(\gamma_{1i}-\pi_{1i})\gamma_{2j}Z_{2j}+\sum\limits_{i,j}a_{ij}\mu_{1i}\pi_{1i}\gamma_{2j}Z_{2j}\\ &&+\sum\limits_{i,j}a_{ij}\mu_{2j}\gamma_{1i}Z_{1i}(\gamma_{2j}-\pi_{2j})+\sum\limits_{i,j}a_{ij}\mu_{2j}\pi_{2j}\gamma_{1i}Z_{1i}\\ &&+\sum\limits_{i,j}a_{ij}\mu_{1i}\mu_{2j}(\gamma_{1i}-\pi_{1i})(\gamma_{2j}-\pi_{2j})+\sum\limits_{i,j}a_{ij}\mu_{1i}\mu_{2j}\pi_{2j}(\gamma_{1i}-\pi_{1i})+\sum\limits_{i,j}a_{ij}\mu_{1i}\pi_{1i}\mu_{2j}(\gamma_{2j}-\pi_{2j})\\ &\equiv&S_{1,1}\\ &&+S_{2,1}+S_{2,2}\\ &&+S_{3,1}+S_{3,2}\\ &&+S_{4,1}+S_{4,2}+S_{4,3}\\ \end{array}

We will use Theorem 1 for the bilinear forms $S_{1,1},S_{2,1},S_{3,1},S_{4,1}$ , and Hoeffding’s inequality, stated below, for the linear combinations $S_{2,2},S_{3,2},S_{4,2},S_{4,3}$ . For $S_{2,1}=\sum\limits_{i,j}a_{ij}\mu_{1i}(\gamma_{1i}-\pi_{1i})\gamma_{2j}Z_{2j}$ , we treat $\gamma_{1i}-\pi_{1i}$ as the centered sub-Gaussian variables with $\psi_{2}$ -norm at most $1$ in which the $n$ sub-Gaussian are completely observed. Applying the union bound and combining the results of $S_{i,j}$ , we can derive the bound for ${\rm P}\Big{[}\big{|}S\big{|}>t\Big{]}$ .

Theorem 5 (Hoeffding’s inequality).

Let $Z_{i}$ , $i=1,\ldots,n$ , be mean-zero independent sub-Gaussian variables satisfying $\max\limits_{1\leq i\leq n}||Z_{i}||_{\psi_{2}}\leq K$ , and $\gamma_{i}$ , $i=1,\ldots,n$ , be independent Bernoulli variables. Also, suppose $Z_{i}$ and $\gamma_{i}$ are independent for all $i$ . Let $\boldsymbol{\alpha}=(\alpha_{i},1\leq i\leq n)^{\top}\in\mathbb{R}^{n}$ be a given coefficient vector. Then, we have

\mathbb{E}\exp\left(\lambda\sum_{i=1}^{n}\alpha_{i}\gamma_{i}Z_{i}\right)\leq\exp\left(C\lambda^{2}K^{2}\sum_{i=1}^{n}\alpha_{i}^{2}\right),\quad\lambda>0,

for some numerical constant $C>0$ . Moreover, we have

{\rm P}\bigg{[}\big{|}\sum\limits_{i=1}^{n}\alpha_{i}\gamma_{i}Z_{i}\big{|}>t\bigg{]}\leq 2\exp\left\{-\dfrac{ct^{2}}{K^{2}||\boldsymbol{\alpha}||_{2}^{2}}\right\},\quad t>0,

for some numerical constant $c>0$ .