The linear inverse problem is given by 𝐘 = 𝐁𝐗 * + 𝐖, where 𝐁 ∈ R 𝑛×𝑑 is the measurement matrix, 𝐘 ∈ R 𝑛×𝑚 represents the linear measurements of the unknown 𝐗 * ∈ R 𝑑×𝑚 , and 𝐖 𝑖𝑗 ∼ (0, 𝜎 2 ) denotes i.i.d. Gaussian noise. In unlabeled sensing, the measurements 𝐘 are scrambled. Specifically,
where 𝐏 * ∈ R 𝑛×𝑛 is the unknown permutation matrix. Given 𝐘 and 𝐁, the objective is to estimate 𝐗 * . In this manuscript, 𝐗 * and 𝐏 * denote the underlying unknown parameters of the unlabeled sensing problem.
If there is no noise, i.e., 𝐖 = 𝟎, we refer to (1) as the exact unlabeled sensing problem. Since estimating 𝐗 * and 𝐏 * for a generic permutation is challenging, several works [1][2][3][4][5][6] assume a 𝑘-sparse or partially shuffled permutation model, Fig. 1. A permutation matrix 𝐏 * 𝑘 is 𝑘-sparse if it has 𝑘 offdiagonal elements, i.e., ⟨𝐈, 𝐏 * 𝑘 ⟩ = 𝑛 -𝑘, where ⟨⋅, ⋅⟩ denotes the trace inner product. The 𝑟-local, or block diagonal, permutation model was first proposed in [7] and later considered in [8,9]. A permutation matrix 𝐏 * 𝑟 is 𝑟-local with 𝑠 blocks if 𝐏 * 𝑟 = blkdiag(𝐏 1 , … , 𝐏 𝑠 ), where blkdiag(⋅) denotes a block diagonal matrix. Fig. 1 illustrates the 𝑘-sparse and 𝑟-local permutation models. These two models are compared by information-theoretic inachievability in [9].
In [9], an alternating minimization algorithm was proposed for the unlabeled sensing problem. The algorithm estimates 𝐏 * and 𝐗 * by minimizing the following optimization program:
( P, X) = argmin 𝐏∈𝛱 𝑛 ,𝐗 𝐹 (𝐏, 𝐗) = ‖𝐘 -𝐏𝐁𝐗‖ 2 𝐹 , (2) where 𝛱 𝑛 denotes the set of 𝑛 × 𝑛 permutations and ‖ ⋅ ‖ 2 𝐹 denotes the squared Frobenius norm of a matrix, which is the sum of the squares of its entries. Given that (2) is a non-convex optimization problem, a crucial part of the alternating minimization algorithm in [9] is the initialization. When the unknown permutation 𝐏 * 𝑟 is 𝑟-local, [9] proposed the collapsed initialization (3). Let 𝐏 𝑖 denote each 𝑟 𝑖 × 𝑟 𝑖 block of 𝐏 * 𝑟 = blkdiag(𝐏 1 , … , 𝐏 𝑠 ). Let 𝐁 𝑖 ∈ R 𝑟 𝑖 ×𝑑 denote blocks of the measurement matrix 𝐁 = [𝐁 1 ; ⋯ ; 𝐁 𝑠 ], where ; denotes vertical concatenation. Then, each block of measurements 𝐘 𝑖 is expressed as 𝐘 𝑖 = 𝐏 𝑖 𝐁 𝑖 𝐗 * , 𝑖 ∈ [𝑠]. The shuffled measurements in each block [𝐏 𝑖 𝐁 𝑖 ∶ 𝐘 𝑖 ] are summed to extract 𝑠 unshuffled measurements [𝟏 ⊤ 𝑟 𝑖
where 𝟏 𝑟 𝑖 ∈ R 𝑟 𝑖 denotes the vector whose entries are all 1. These 𝑠 measurements are then represented compactly in the collapsed linear system of equations B𝐗 * = Ỹ, where B ∈ R 𝑠×𝑑 and Ỹ ∈ R 𝑠×𝑚 . The initialization for 𝑟-local
𝑟 is the minimum-norm solution to the collapsed system and is given by:
https://doi.org/10.1016/j.sigpro.2025.109927 Received 23 July 2024; Received in revised form 11 December 2024; Accepted 27 January 2025 Left. Sparse (or partially shuffled) permutation considered in [1][2][3][4], with number of shuffles 𝑘 = 10. Right. The 𝑟-local permutation structure considered in [7,9,10], with block size 𝑟 = 10. In this paper, we propose a general algorithm for both permutation models.
where X(0
For the 𝑘-sparse problem, in this manuscript, we propose the following initialization:
Note that the above initialization corresponds to the least square solution of (1) when 𝐏 * is the identity matrix. The initialization is motivated by the observation that, since 𝐏 * 𝑘 has 𝑘 off-diagonal elements, the identity matrix serves as a simple approximation of 𝐏 * 𝑘 . The main goal of this manuscript is to characterize the effectiveness of the initializations 𝐗 (0) 𝑟 and 𝐗 (0) 𝑘 by upper bounding the initialization error for the exact unlabeled sensing problem.
This paper presents a theoretical analysis of the initialization for the exact unlabeled sensing problem under two permutation models. The key contributions of this manuscript are summarized as follows:
1. Assuming 𝐁 is Gaussian, we show that the relative error
critically depends on the parameter 𝑑 -𝑠. In particular, Theorem 4.1 establishes that the probability that the relative error exceeds √ 1 -𝑠∕𝑑 decays at a sub-Gaussian rate.
2. Assuming 𝐱 * is sub-Gaussian, we focus on the error
, which also depends critically on the parameter 𝑑 -𝑠. In particular, Theorem 4.2 in this work shows that the probability that the error exceeds 𝑂 ( √ 𝑑 -𝑠 ) also decays at a sub-Gaussian rate. 3. For the 𝑘-sparse problem, we analyze the relative error
. Assuming 𝐁 is Gaussian and sufficiently tall (see Definition 4.4), we show that the relative error depends on the parameter 𝑘∕𝑛. Specifically, depending on the problem parameters, Theorem 4.5 shows that the probability of the error exceeding 𝑂(𝑘∕𝑛) decays at a rate ranging from inverse linear to exponential. 4. We apply the alternating minimization algorithm in [9] with the initializations (3) and ( 4). Experiments results on real datasets show superior performance compared to baseline methods.
Outline. The outline of the paper is as follows. Section 2 discusses related work. Section 3 presents the alternating minimization algorithms for the 𝑟-local and 𝑘-sparse permutation models. Section 4 covers the initialization analysis. Experimental results on real datasets showing superior performance compared to baseline methods are given in Section 5. We conclude in Section 6.
Table 1
Table summarizing the frequently used notation used in this paper. See Fig. 1 for definition and examples of 𝑟-local 𝐏 * 𝑟 and 𝑘-sparse 𝐏 * 𝑘 permutations. The superscript * in 𝐏 * , 𝐏 * 𝑟 , 𝐏 * 𝑘 , 𝐗 * denotes that these are unknown problem parameters which are estimated by the proposed algorithm. Symbol Description 𝐁 ∈ R 𝑛×𝑑 Measurement matrix 𝐗 * ∈ R 𝑑×𝑚 Unknown matrix 𝐱 * 𝑗 ∈ R 𝑑 𝑗-th column of 𝐗 * , ∀𝑗 ∈ [𝑚] 𝐏 * ∈ {0, 1} 𝑛×𝑛 Unknown permutation matrix 𝐏 * 𝑘 ∈ R 𝑛×𝑛 Unknown 𝑘-sparse permutation matrix 𝑘 Number of shuffles 𝑘 s.t. ⟨𝐏 * 𝑘 , 𝐈⟩ = 𝑛 -𝑘 𝐏 * 𝑟 ∈ {0, 1} 𝑛×𝑛 Unknown 𝑟-local permutation matrix 𝑟 𝑖 ∀ 𝑖 ∈ [𝑠] Block sizes of 𝐏 * 𝑟 s.t. 𝐏 * 𝑟 = blk diag (𝐏 * 1 , … , 𝐏 * 𝑠 ) 𝛱 𝑟 𝑖 The set of 𝑟 𝑖 × 𝑟 𝑖 permutation matrices. 𝑠 Number of blocks in 𝑟-local 𝐏 * 𝑟 [𝑠] Integers {1, … , 𝑠}
Boldface lowercase letters denote column vectors and boldface uppercase letters denote matrices. The transpose of a vector 𝐱 is noted by 𝐱 ⊤ . For a matrix 𝐗, its transpose is denoted by 𝐗 ⊤ . 𝐀 † denotes the Moore-Penrose pseudoinverse of 𝐀. Tr (⋅) denotes the trace of a matrix. Given two matrices 𝐀 and 𝐁, ⟨𝐀 , 𝐁⟩ denotes the trace inner product, defined as Tr (𝐀 ⊤ 𝐁). Given a vector 𝐱, the 𝑖-th entry of 𝐱 is denoted by 𝑥 𝑖 . [𝐀; 𝐁] denotes the vertical concatenation of matrices 𝐀, 𝐁. ‖ ⋅ ‖ 𝑝 denotes the vector 𝑝-norm. ‖ ⋅ ‖ 𝐹 and ‖ ⋅ ‖ denote the Frobenius norm and the operator norm of a matrix, respectively. 𝟏 𝑛 ∈ R 𝑛 denotes the vector of all ones. 𝐈 denotes the identity matrix. R denotes the set of real numbers. 𝛱 𝑛 = {𝐙 ∶ 𝐙 ∈ {0, 1} 𝑛×𝑛 , 𝐙𝟏 𝑛 = 𝟏 𝑛 , 𝐙 ⊤ 𝟏 𝑛 = 𝟏 𝑛 } denotes the set of 𝑛 × 𝑛 permutation matrices. E(⋅) denotes the expectation of a random variable in consideration. 𝐶 , 𝐶 ′ denote absolute constants > 1, and 𝑐 , 𝑐 ′ denote absolute constants ≤ 1. argmin denotes the set of minima of an objective function in consideration. exp(⋅) denotes the exponential function. We also summarize the frequently used notation in Table 1.
This section provides a concise overview of related work in unlabeled sensing theory and algorithms, inference problems involving unlabeled data, and applications of unlabeled sensing.
We note that the problem in (1) is referred to as the single-view (multi-view) unlabeled sensing problem for 𝑚 = 1 (𝑚 > 1). For the single-view problem, we represent the problem as 𝐲 = 𝐏 * 𝐁𝐱 * + 𝐰 where 𝐲, 𝐱 * , and 𝐰 denote the variables analogous to 𝐘, 𝐗 * and 𝐖 respectively. The work in [11] formulated the single-view unlabeled sensing problem and established that 𝑛 = 2𝑑 measurements are both necessary and sufficient for the recovery of 𝐱 * . Subsequent works [12,13] generalized this result and developed an information-theoretic inachievability analysis.
For single-view unlabeled sensing, algorithms based on branch and bound and expectation maximization are proposed in [14][15][16], which are suitable for small problem sizes. A stochastic alternating minimization approach is introduced in [17]. For multi-view unlabeled sensing, the Levsort subspace matching algorithm was proposed in [18], and a subspace clustering approach was presented in [4]. The works [1,2,19] propose methods based on bi-convex optimization, robust regression, and spectral initialization, respectively.
The problem of unlabeled sensing has been studied in the context of estimation and detection from unlabeled samples. In [20], the authors consider the problem of signal amplitude estimation and detection from unlabeled quantized binary samples. In the setup they consider, the ordering of the time indexes is unknown. The work in [20] proposes a maximum likelihood estimator to estimate the underlying signal amplitude and permutation, under certain structural assumptions on the quantization and signal profile. Furthermore, an alternating maximization algorithm is studied for the general estimation and detection problem.
In [21], the authors study signal recovery from unlabeled samples by considering a special case of unlabeled sensing, referred to as the unlabeled ordered sampling problem. In this problem, instead of estimating an unknown permutation matrix, the goal is to estimate an unknown selection matrix that preserves the order of the measurements. [21] links this problem to compressive sensing and also proposes an alternating minimization algorithm for solving it. [22] also addresses the unlabeled ordered sampling problem, but with a focus on signal detection, rather than signal recovery as in [21]. A dynamic programming algorithm is proposed therein to estimate the selection matrix.
[23] considers the problem of signal detection, where the observations are independently drawn from a finite alphabet. Focusing on the inference problem, the work in [23] characterizes the information available in the unordered samples by studying a binary hypothesis test. Additionally, it provides a computationally efficient detector to address the detection problem. [24] further studies the signal detection problem in the case of binary unlabeled observations. In the low-detectability regime, [24] gives an analytical characterization of various statistical inference quantities, such as Chernoff information. In [25], the signal detection framework of unlabeled sensing is applied to decision-making in sensor networks, where sensors report measurements to a central node, but the measurements are imprecise or lack labels. [26] considers a related problem, where sensors send their measurements to a central fusion center. Due to an anomaly in one of the sensors, the measurements at the fusion center are assumed to be unordered, and [26] explores the problem of detecting the anomaly with minimal detection delay.
The linked linear regression problem [27][28][29], also called regression analysis with linked data, is to fit a linear model by minimizing over 𝐱 and 𝐏 the objective ‖𝐲 -𝐏𝐁𝐱‖ 2 , where 𝐏 is an 𝑟-local permutation and 𝐱 is the regression vector. For example, consider the statistical model where the weight depends linearly on age and height. Let 𝐲 ∈ R 𝑛 contain the weights of 𝑛 = 10 individuals, 5 males, 5 females. Let 𝐁 ∈ R 𝑛×𝑑 , with 𝑑 = 2 contain the age and height values. Each record (weight, age, height) is known to belong to a male or a female individual. Letting the first (second) block of 5 rows of 𝐲, 𝐁 contain the records of male (female) individuals, the unknown 𝑟-local permutation, 𝑟 = 5, assigns each weight value (row of 𝐲) to its corresponding age, height (row of 𝐁). Detailed references describing record linkage with blocking are [30], Section 1.1 and [31], Section 2; also, see Fig. 2. Experimental results on real datasets are given in Section 5. Other applications of unlabeled sensing include the pose and correspondence problem [13] and 1-d unassigned distance geometry [9]. In [32], the authors apply the unlabeled sensing framework to the sensor network localization problem. We also note a recent work that extends the unlabeled sensing theory to multi-channel signals, leading to highly structured unlabeled sensing problems [33]. In [33], beyond theoretical analysis of the structured problem, the model is applied to a real-world application in calcium imaging.
This section provides definition of sub-Gaussian and sub-exponential random vectors, and also summarizes the key concentration inequalities used throughout this paper. Let 𝐙 * ∈ R 𝑑×𝑚 be a matrix such that the columns 𝐳 * 𝑗 are independent, zero-mean sub-Gaussian random vectors, with sub-Gaussian norm 𝐾. It follows then that
A random variable 𝑋 is sub-exponential with sub-exponential norm
where 𝐶 ≥ 1, 𝑡 ≥ 0 and 𝐾 𝑑-𝑠 ≥ 0.
Theorem 2.1 (Hanswon-Wright Inequality, Theorem 2.1 in [35]). Let Σ = 𝐀 ⊺ 𝐀 be a positive semi-definite matrix. Let 𝐱 = (𝑥 1 , … , 𝑥 𝑑 ) be a zero-mean sub-Gaussian random vector, i.e., for
Lemma 2.2 (Johnson-Lindenstrauss Lemma, Lemma 5.3.2 in [36]). Let 𝑃 be a projection in R 𝑝 onto a uniformly distributed random 𝑞-dimensional subspace. Let 𝐳 ∈ R 𝑝 be a fixed point and 𝑡 > 0. Then, with probability at least 1 -2 exp(-𝑐 𝑡 2 𝑞),
Theorem 2.3 (Hoeffding's Inequality, Theorem 2.6.2 in [36]). Let 𝑋 1 , … , 𝑋 𝑚 be independent, sub-Gaussian random variables. Then, for every 𝑡 ≥ 0,
where ‖𝑋 𝑖 ‖ 𝜑 2 denotes the sub-Gaussian norm of 𝑋 𝑖 and 𝑐 ′ is an absolute constant.
Lemma 2.4 (Tail Inequality for 𝜒 2 𝐷 Distributed Random Variables, Lemma 1 in [37]). Let 𝑍 𝐷 be a 𝜒 2 statistic with 𝐷 degrees of freedom. For any positive 𝑡,
In this section, we briefly summarize the alternating minimization algorithm first proposed in [9]. To estimate the 𝐏 * and 𝐗 * in the optimization objective (1), the alternating minimization (AltMin) updates for (2) are
where 𝐏 (𝑡) and 𝐗 (𝑡) denote the estimate of 𝐏 * and 𝐗 * at the t-th iteration respectively. The initialization to ( 12), ( 13), a linear assignment problem and a least squares problem, has to be chosen carefully because ( 2) is a non-convex optimization problem. After initialization, we alternate ( 12), ( 13) and use the relative change in the objective value 𝐹 (𝐗, 𝐏) as the stopping criterion. For 𝑟-local 𝐏 * 𝑟 , the 𝐏-update (12) decouples along the blocks of 𝐏 * 𝑟 , see Algorithm 1. For 𝑘-sparse 𝐏 * 𝑘 , see Algorithm 2.
Fig. 2. Record linkage with Blocking [27,34] assigns records in File A to records in File B, up to blocks. For example, records matching on identifiers 'city', 'occupation' etc. (left) are assigned to the same block (right). Linked linear regression [28] fits a regression model on such block-permuted data. See Section 5 for results on real datasets. The figure on the left is adapted from Figure 1 in [34] and the figure on the right is adapted from Figure 1 in [27].
Algorithm 1 AltMin for 𝑟-local 𝐏 * 𝑟 Input: Convergence threshold 𝜖, 𝐘, 𝐁 and block sizes 𝑟 1 , 𝑟 2 , .., 𝑟 𝑠 of 𝐏 * 1 , 𝐏 * 2 , ..., 𝐏 * 𝑠 respectively. 1: X ← collapsed initialization in (3) 2: Ŷ ← 𝐁 X 3: while relative change > 𝜖 do 4: for 𝑖 ∈ 1 ⋯ 𝑠 do //𝑠 is the number of blocks 5: P𝑖 ← argmin 𝐏 𝑖 ∈𝛱 𝑟 𝑖 -⟨𝐘 𝑖 Ŷ⊤ 𝑖 , 𝐏 𝑖 ⟩ 6: P ← blkdiag( P1 , ⋯ , P𝑠 ) 7: X ← 𝐁 † P⊤ 𝐘 8: Ŷ ← 𝐁 X 9: Return P, X Algorithm 2 AltMin for 𝑘-sparse 𝐏 * 𝑘 Input: convergence threshold 𝜖, 𝐘, 𝐁 1: Ŷ ← 𝐘 2: while relative change > 𝜖 do 3: P ← argmin 𝐏∈𝛱 𝑛 -⟨𝐘 Ŷ⊤ , 𝐏⟩ 4:
Ŷ ← 𝐁 X 6: Return P, X
In this section, we derive upper bounds for the initialization error of the alternating minimization algorithm applied to solve the exact unlabeled sensing problem ((1) with 𝐖 = 𝟎), using the initializations in (3) and (4).
We consider the initialization in (3). Let B = Ũ𝐒 Ṽ⊤ denote the compact singular value decomposition of B. Using this in (3), we obtain B † = Ṽ𝐒 -1 Ũ⊤ , where B = Ũ𝐒 Ṽ⊤ . Using this substitution, it can be shown that the initialization error ‖𝐗 * -X(0) 𝑟 ‖ 𝐹 is the projection of 𝐗 * onto the orthogonal complement of the row space of B, as follows:
Recall that the size of B is 𝑠 × 𝑑. If 𝑠 ≥ 𝑑, and assuming that B has rank 𝑑, it can be verified that X(0) 𝑟 = 𝐗 * . For instance, for 𝑠 ≥ 𝑑, if the entries of 𝐁 are drawn independently from a continuous distribution, the rank of B is 𝑑 with high probability. Given that, in this section, we upper bound the error in initialization for the under-determined 𝑠 < 𝑑 case. Similar analysis for under-determined systems have been studied in the sketching literature [38], but these results are not applicable here as our sub-sampling strategy via the collapsing initialization is deterministic.
The first key theoretical result is Theorem 4.1. This theorem considers a Gaussian matrix 𝐁 ∈ R 𝑛×𝑑 , with 𝐗 * ∈ R 𝑑×𝑚 assumed to be a fixed but unknown matrix. To upper bound the relative error
we apply the Johnson-Lindenstrauss lemma. The key insight is that the relative error depends on the parameter 𝑑 -𝑠. Specifically, the probability that the relative error exceeds √ 1 -𝑠 𝑑 decays at a sub-Gaussian rate. This implies that the initialization in (3) improves as the number of blocks 𝑠 in the 𝑟-local model increases. 𝑟 be as defined in (3). Consider the exact unlabeled sensing problem. Assuming Gaussian 𝐁 ∈ R 𝑛×𝑑 and block-diagonal 𝐏 * 𝑟 with 𝑠 blocks, for log 𝑚 ≤ 𝑐 𝑡 2 (𝑑 -𝑠), 𝑠 < 𝑑, and 𝑡 ≥ 0, Pr
Proof. The error ( 14) is the column-wise projection of 𝐗 * ∈ R 𝑑×𝑚 onto a (𝑑 -𝑠)-dimensional uniformly random subspace (14), which can be bounded by the JL Lemma (8). □
In Theorem 4.1, we considered a Gaussian matrix 𝐁 and a fixed, unknown matrix 𝐗 * . In the following, we fix the measurement matrix 𝐁 and introduce randomness in 𝐗. Before discussing our main result for this setting, we provide motivation for considering this scenario. One motivating application is the unassigned distance geometry problem (uDGP) [39][40][41], which involves recovering the configuration of points from a list of distances. We focus on the 1-dimensional uDGP, which can be framed as a structured unlabeled sensing problem, as described in [9], where 𝐱 * represents the sought-after positions of the points on the 1D line. In this case, the matrix 𝐁 is deterministic (see Eq. ( 3) of [9]). For this setting, without any assumptions on the underlying points 𝐱 * , and assuming that the underlying permutation is 𝑟-local, the collapsed initialization may be suboptimal. To illustrate this, recall the error in the collapsed estimate: ‖(𝐈 -Ṽ Ṽ⊤ )𝐱 * ‖ 𝐹 . For any 𝐱 * orthogonal to the subspace spanned by Ṽ, the collapsed estimate will be the zero vector. This shows that the following upper bound on the initialization error ‖ x(0) 𝑟 -𝐱 * ‖ ≤ ‖𝐈 -Ṽ Ṽ⊤ ‖‖𝐱 * ‖ = ‖𝐱 * ‖, holds with equality for 𝐱 * orthogonal to the subspace spanned by the columns of Ṽ. The implication of the above discussion is that some structural assumption on the underlying 𝐱 * is necessary to obtain suitable bounds on the initialization error.
We next consider the squared error in the initialization ‖𝐀𝐱 * ‖ 2 where 𝐀 = 𝐈 -Ṽ Ṽ⊤ . Note that 𝐀 is an orthogonal projection matrix onto the orthogonal complement of the subspace spanned by Ṽ. Given that Ṽ ∈ R 𝑑×𝑠 and we are considering the underdetermined regime 𝑠 < 𝑑, the rank of 𝐀 is 𝑑 -𝑠. Therefore, the analysis of the error is equivalent to analyzing the quadratic form ‖𝐀𝐱 * ‖ 2 . One common approach to control the norm of a random vector is to assume that 𝐱 * is a random vector with independent, sub-Gaussian coordinates (see, for example, Theorem 3.1.1 in [36]). However, we note that this assumption does not hold after applying the projection operator 𝐀. A more suitable tool in this case is the Hanson-Wright inequality. In fact, the quadratic form in our case is special due to the fact that 𝐀 is a projection operator, allowing us to leverage its linear algebraic properties. To apply the Hanson-Wright inequality, it suffices to assume that 𝐱 * is sub-Gaussian with mean zero, without requiring the independence of its coordinates. With this assumption, we study the error
4.1 shows that this error depends critically on the parameter 𝑑-𝑠. In particular, Theorem 4.2 demonstrates that the probability of the error exceeding 𝑂 ( √ 𝑑 -𝑠 ) decays at a sub-Gaussian rate. Theorem 4.2. Let 𝐗 * ∈ R 𝑑×𝑚 be the unknown matrix such that the columns 𝐱 * 𝑗 are independent, zero-mean sub-Gaussian random vectors, with sub-Gaussian norm 𝐾. Let X(0) 𝑟 = [ x(0) 1 | ⋯ | x(0) 𝑚 ] be as defined in (3). Consider the exact unlabeled sensing problem. For any fixed measurement matrix 𝐁 and block-diagonal 𝐏 * 𝑟 with 𝑠 blocks, 𝑠 < 𝑑 and 𝑡 ≥ 0,
where
The strategy for the proof of Theorem 4.2 is to first use the Hanson-Wright inequality to control the quadratic form ‖𝐀𝐱 * ‖ 2 . Once that is established, we use Lemma 4.3 to define a modified random variable that is sub-exponential. The rest of the technical proof uses the Hoeffding bound on the modified random variable and standard techniques to relate the concentration inequality of ‖𝐳‖ 2 2 to ‖𝐳‖ 2 [36]. For the complete proof, see Appendix A.
𝑟 be as in (3), 𝐱 * ∈ R 𝑑 be a zero-mean sub-Gaussian random vector, 𝐁 be a fixed measurement matrix, and 𝐏 * 𝑟 be a fixed blockdiagonal permutation with 𝑠 blocks. Let 𝐱 err denote the random variable
1. For 𝑠 < 𝑑 and 𝑡 ≥ 0,
where
The random variable xerr is a sub-exponential random variable with sub-exponential norm ‖x err -𝐶 𝑑-𝑠 ‖ 𝜑 1 = 𝐶 𝐾 𝑑-𝑠 . In addition, xerr ≥ 𝐱 err , and
Proof. We prove each part separately.
1. To derive (17), set 𝐀 = (𝐈 -Ṽ Ṽ⊤ ) in (7), where Ṽ ∈ R 𝑠×𝑑 denotes the basis for the row space of the collapsed matrix
Here, 𝐀 ∈ R 𝑑×𝑑 is a (𝑑 -𝑠)-dimensional projection matrix, and Σ = 𝐀 ⊤ 𝐀 = 𝐀𝐀 = 𝐀. Moreover, Tr (Σ) = Tr (Σ 2 ) = Tr (𝐀) = -𝑠 and the operator norm satisfies ‖Σ‖ = 1. We now apply Hanson-Wright inequality in (7), we obtain:
The next step will be to upper bound the term 2 √ (𝑑 -𝑠)𝑡. To do that, we rely on following inequality: 2 √ 𝑏𝑡 ≤ 2𝑡 + 1 2 𝑏, for 𝑏 ≥ 1 and 𝑡 ≥ 0. To verify this, note that this inequality is equivalent to checking whether 4𝑡 ≤ 4𝑡 2 𝑎 + 2𝑡𝑎 + 𝑎∕4 holds.
• For 𝑎 = 1, this reduces to 0 ≤ (2𝑡 -1∕2) 2 , which holds for all 𝑡 ≥ 0. • For 𝑎 ≥ 2, we check 4𝑡 ≤ 4𝑡 2 𝑎 + 2𝑡𝑎 + 1∕4 which is true for 𝑡 ≥ 0.
Using this inequality with 𝑏 = 𝑑 -𝑠, we obtain 2 √ (𝑑 -𝑠)𝑡 ≤ 2𝑡(𝑑 -𝑠) + 1 2 (𝑑 -𝑠). Substituting this bound in (21), we simplify the inequality as follows:
The above inequality is equivalent to:
A change of variable, 𝑢 = 2𝐾 2 𝑡((𝑑 -𝑠) + 1), and minor algebraic manipulation yields the final desired result.
Remark: To derive Theorem 4.2, we first show that 𝐱 err is a subexponential random variable. It cannot be concluded from ( 17) that the shifted random variable 𝐱 err -𝐶 𝑑-𝑠 is sub-exponential because the lower tail probability Pr [𝐱 err -𝐶 𝑑-𝑠 ≤ 0] is not bounded. Since our goal is to upper bound the probability that the error exceeds a certain value, we define the non-negative subexponential random variable xerr -𝐶 𝑑-𝑠 , which upper bounds the lower tail as {𝐱 err ≤ 𝐶 𝑑-𝑠 } = 𝐶 𝑑-𝑠 , see (18). We will use the definition of xerr given in ( 18), ( 19). 2. Conditioning on the two complementary events in ( 18), (19), by the law of total probability, ∀ 𝑡 > 0, we have
= Pr [𝐱 err -𝐶 𝑑-𝑠 ≥ 𝑡 | 𝐱 err > 𝐶 𝑑-𝑠 ] Pr [𝐱 err > 𝐶 𝑑-𝑠 ] (25) = Pr [𝐱 err -𝐶 𝑑-𝑠 ≥ 𝑡], (26) (24) follows from (18); specifically, Pr [x err -𝐶 𝑑-𝑠 ≥ 𝑡 | 𝐱 err -𝐶 𝑑-𝑠 ] = 0 ∀ 𝑡 > 0. (25) follows from (19). (26) follows from noting that the event {𝐱 * | 𝐱 err ≥ 𝐶 𝑑-𝑠 + 𝑡} is a subset of {𝐱 * | 𝐱 err > 𝐶 𝑑-𝑠 }. From (26), for 𝑡 > 0 Pr [𝐱 err -𝐶 𝑑-𝑠 ≥ 𝑡] = Pr [|x err -𝐶 𝑑-𝑠 | ≥ 𝑡]. (27) (27) follows from noting that xerr -𝐶 𝑑-𝑠 = |x err -𝐶 𝑑-𝑠 |. Substituting (17) in (27), for 𝑡 > 0, Pr [|x err -𝐶 𝑑-𝑠 | ≥ 𝑡] ≤ exp ( -𝑐 𝑡 𝐾 𝑑-𝑠
Since exp(0) = 1, (28) holds for 𝑡 ≥ 0. In ( 28), we have verified the definition in ( 6) for xerr -𝐶 𝑑-𝑠 with 𝐾 𝑑-𝑠 = 2𝐾 2 ( √ 𝑑 -𝑠 + 1) = 𝐶 𝐾 2 √ 𝑑 -𝑠. By proposition 2.7.1 (a), (d), definition 2.7.5 in [36] and ( 28), the sub-exponential norm is ‖x err -𝐶 𝑑-𝑠 ‖ 𝜑 1 = 𝐶 𝐾 𝑑-𝑠 . □
In this section, we analyze the initialization step when the underlying permutation is 𝑘-sparse. Specifically, we study the relative error
. Assuming that 𝐁 is a ''tall'' Gaussian matrix, Theorem 4.5
establishes that this error depends on the parameter 𝑘∕𝑛, which represents the proportion of shuffled entries relative to the total. More precisely, Theorem 4.5 shows that the probability of the error exceeding 𝑂(𝑘∕𝑛) decays at a rate ranging from inverse linear to exponential. Central to Theorem 4.5 is the assumption that the measurement matrix 𝐁 is a ''tall'' Gaussian matrix, which we define precisely below.
). A Gaussian matrix 𝐁 ∈ R 𝑛×𝑑 is considered ''tall'' if the aspect ratio 𝜆 = 𝑑∕𝑛 satisfies 𝜆 < 𝜆 0 for some sufficiently small constant 𝜆 0 > 0. In that case, we have Pr [𝜎 min (𝐁) ≥ 𝑐 √ 𝑛] ≤ 𝑒 -𝑐 𝑛 , where 𝑐 is an absolute constant that depends on 𝐁. , where 𝐶 = 4(1 + √ 3), and 1 ≤ 𝑘 ≤ 𝑛, we have
For the proof of ( 29), see Appendix A in the Appendix.
In the above theorem, we consider the error bound at the lower bounds and upper bounds of 𝜖.
, it can be shown that the failure probability is 8
, the failure probability is exp(-𝑐 𝑛). (29) gives a (< 1) relative error bound when 𝑘 grows slowly with 𝑛, for example 𝑘 = 𝑛 𝛽 , 𝛽 < 1.
We compare the proposed AltMin algorithm to several benchmark methods on both synthetic (Fig. 3) and real (Table 3) datasets. We implemented all algorithms in MATLAB. The linear assignment problem to recover the permutation estimate 𝐏 is solved by using the MATLAB matchpairs function. The least squares estimate (line 7 of Algorithm 1 and line 4 of Algorithm 2) is solved by computing the Moore-Penrose pseudoinverse using the MATLAB function pinv. The convex optimization problem proposed in [2] is solving using CVX [43,44]. The specific solver we used is SeDuMi. We also use CVX to project onto the set of doubly stochastic matrices for the benchmark method 'ds+'. Baselines. We compare against six baseline methods, '𝓁 2 -regularized' [2], 'ds+' [2], 'Spectral' [19], 'Biconvex' [1], 'RLUS' [7], and 'Stochastic AltMin' [17]. The '𝓁 2 -regularized' method considers the 𝑘-sparse permutation model and imposes a row-wise group sparse penalty (𝓁 2 -norm regularization) on 𝐌, where Ŷ = 𝐁 X + 𝐌. Other methods are discussed in the following paragraphs. For the proposed algorithms, ( 12), ( 13) are alternated until the change in the objective value is less than 1 percent.
To simulate 𝐘 = 𝐏 * 𝐁 𝑛×𝑑 𝐗 * 𝑑×𝑚 + 𝐖, we generate data by drawing the entries of matrices 𝐁, 𝐗 * and 𝐖 from the normal distribution. Subsequently, 𝐖 is scaled by 𝜎 to set SNR ≜ ‖𝐗 * ‖ 2 𝐹 ∕(𝑚𝜎 2 ) = 100. The permutation matrices 𝐏 * 𝑟 and 𝐏 * 𝑘 are sampled uniformly from the set of 𝑟-local and 𝑘-sparse permutations, respectively. The results are averaged over 15 Monte-Carlo runs.
To adapt the 'Spectral', '𝓁 2regularized', and 'Biconvex' methods to the 𝑟-local model, we add a constraint enforcing the permutation estimates to be block-diagonal.'
• The Spectral method in [3] studies the unlabeled sensing model 𝐘 = 𝐏 𝑖 𝐁𝐗 + 𝐖, where 𝐖 is an additive Gaussian noise, with no assumption made on the structure of the permutation, i.e., the underlying permutation is generic. The proposed algorithm in [3] Table 3 𝑟 max denotes the largest block size of 𝐏 * 𝑟 and 𝑠 denotes the number of blocks. For a description of the datasets, see Table 2. The methods are compared by the coefficient of determination 𝑅 2 and relative error (𝑅 2 , relative error). The relative error is ‖𝐗 * -X‖ 𝐹 ∕‖𝐗 * ‖ 𝐹 , where 𝐗 * = 𝐁 † 𝐘 * is the 'oracle' regression matrix given unpermuted data 𝐘 * . The coefficient 𝑅 2 ( X) = 1 -(‖𝐘 * -𝐁 X‖ 𝐹 ∕‖𝐘 * ‖ 𝐹 ) measures the goodness of fit for the unpermuted data. The 'naive' estimate from permuted data 𝐘 is X = 𝐁 † 𝐘. The coefficient 𝑅 2 is bounded 0 ≤ 𝑅 2 ≤ 1, and a higher value of 𝑅 2 indicates a better fit.
𝑟 max 𝑠 Oracle Naive RLUS 𝓁 2 -reg AltMin-𝑘 AltMin-𝑟 ftp 43 46 (0.88, 0) (0.44, 0.70) (0.62, 0.01) (0.57, 0.59) (0.73, 0.41) (𝟎.𝟖𝟓, 𝟎.𝟏𝟕) ames 311 61 (0.85, 0.32) (0.37, 0.75) (0.38, 0.72) (0.72, 0.89) (0.69, 0.42) (𝟎.𝟕𝟔, 𝟎.𝟑𝟐) scm 424 1137 (0.58, 0) (0.21, 0.79) (0.53, 0.29) (0.46, 0.49) (0.54, 0.29) (𝟎.𝟓𝟓, 𝟎.𝟐𝟓) scs 3553 356 (0.81, 0) (0.01, 0.99) (0.74, 0.24) (0.02, 0.98) (0.73, 0.33) (𝟎.𝟕𝟕, 𝟎.𝟐𝟏) air-qty 95 366 (0.69, 0) (0.29, 0.78) (0.65, 0.20) (0.64, 0.21) (𝟎.𝟔𝟔, 𝟎.𝟏𝟕) (0.65, 0.19) air-qty 46 744 (0.69, 0) (0.10, 0.94) (0.60, 0.34) (0.26, 0.79) (0.57, 0.42) (𝟎.𝟔𝟐, 𝟎.𝟑𝟐)
estimates 𝐏 using
In the case that the underlying permutation is 𝑟-local, first recall that 𝑠 denotes the number of blocks in the 𝑟-local permutation and 𝑟 𝑖 is the size of the 𝑖th block. The unlabeled sensing model reduces to
where 𝐘 𝑖 ∈ R 𝑟 𝑖 ×𝑚 , 𝐁 𝑖 ∈ R 𝑟 𝑖 ×𝑑 and 𝐖 𝑖 ∈ R 𝑟 𝑖 ×𝑚 respectively denote blocks of the matrices 𝐘 ∈ R 𝑛×𝑚 , 𝐁 ∈ R 𝑛×𝑑 and 𝐖 ∈ R 𝑛×𝑚 . Therefore, accounting for the structure of the 𝑟-local permutation, we modify (30) as follows:
which are 𝑠 linear assignment problems over the sets of permutation matrices of size 𝑟 𝑖 . • For the 'Biconvex' algorithm, we modify the 𝐏 1 , 𝐏 2 updates in Algorithm 1 from [1].
where 𝜇 and 𝜌 are ADMM parameters, and 𝐏 𝐁 is the projection matrix onto the column-span of 𝐁. Instead of updating 𝐏 (𝑡+1) 1 = argmin 𝐏∈𝛱 𝑛 ⟨𝐏, 𝐂 1 ⟩, we enforce the block diagonal constraint by updating 𝐏 (𝑡+1) 1,𝑖 = argmin 𝐏 𝑖 ∈𝛱 𝑟 𝑖 ⟨𝐏 𝑖 , 𝐂 1,𝑖 ⟩, where 𝐂 1,𝑖 ∈ R 𝑟×𝑟 is the matrix comprised of entries from the 𝑟 𝑖 rows and 𝑟 𝑖 columns in the 𝑖-the block of 𝐂. We update 𝐏 (𝑡+1) 2 in a similar manner. • The '𝓁 2 -regularized' method considers the 𝑘-sparse permutation model and imposes a row-wise group sparse penalty on 𝐌, where Ŷ = 𝐁 X + 𝐌. Specifically, the proposed estimate is obtained by solving the minimization problem: min 𝐗,𝐌 ‖𝐘 -𝐁𝐗‖ 2 𝐹 + 𝜆 𝑛 ∑ 𝑖=1 ‖𝐌 𝑖 ‖ 2 , where 𝐌 𝑖 denotes the 𝑖th row of 𝐌. We do not modify this minimization over 𝐗 and 𝐌. Instead, We substitute this estimate, denoted by X, and let Ŷ = 𝐁 X. To estimate the underlying permutation, we consider the following minimization problem min 𝐏 ‖𝐘 -𝐏 Ŷ‖ 2 𝐹 . We note the following equalities: argmin 𝐏 ‖𝐘 -𝐏 Ŷ‖ 2 𝐹 = argmin 𝐏 ‖𝐘‖ 2 𝐹 + ‖𝐏 Ŷ‖ 2 𝐹 -2⟨𝐘, 𝐏 Ŷ⟩ = argmax 𝐏 ⟨𝐘, 𝐏 Ŷ⟩ = argmax 𝐏 Tr (𝐘 ⊤ 𝐏 Ŷ) = argmax 𝐏 ⟨𝐘 Ŷ⊤ , 𝐏⟩. In the case that the underlying permutation is 𝑟-local, the above minimization problem decouples, and we obtain a block-wise linear assignment problem: P𝑖 ← argmin 𝐏 𝑖 ∈𝛱 𝑟 𝑖 -⟨𝐘 𝑖 Ŷ⊤ 𝑖 , 𝐏 𝑖 ⟩, where Ŷ𝑖 ∈ R 𝑟 𝑖 ×𝑚 and 𝐘 𝑖 ∈ R 𝑟 𝑖 ×𝑚 denote blocks of the matrices Ŷ and 𝐘, respectively. conclusions of Theorems 4.1, 4.2, and 4.5 as the initialization improves with lower values of 𝑟 and 𝑘.
The proposed AltMin algorithm is also applicable to both models and computationally scalable. For 𝑟 = 125 (Fig. 3(a)) and 𝑘 = 850 (Fig. 3(b)), MATLAB runtime with 16 Gb RAM and 9-th Gen. 4-core processor is less than a second. Fig. 3(c). The entries of the measurement matrix 𝐁 are sampled i.i.d. from the uniform [0, 1] distribution. Compared to the case for Gaussian 𝐁 (Fig. 3(b)), the performance of the 'Spectral' and 'RLUS' methods deteriorates significantly. This is because both algorithms consider quadratic measurements 𝐘𝐘 ⊤ . Specifically, the 'Spectral' method is based on the spectral initialization technique [45], which assumes a Gaussian measurement matrix. In contrast, the performance of the proposed and '𝓁 2 -regularized' methods does not deteriorate. Fig. 3(d). The 'ds+' algorithm [2] considers the convex relaxation of (2) by minimizing the objective over the set of doubly stochastic matrices. Assuming a known upper bound on the number of shuffles 𝑘, 'ds+' constrains ⟨𝐈, 𝐏⟩ ≥ 𝑛 -𝑘. To project onto the set of doubly stochastic matrices, each iteration of 'ds+' minimizes a linear program which greatly increases its run-time. The proposed AltMin algorithm optimizes the same objective, but over the set of permutation matrices. This results in a simpler linear assignment problem, which is why AltMin is faster and outperforms 'ds+'. Fig. 3(e). We compare to the method in [17] which considers the 𝑚 = 1 single-view setup and proposes stochastic alternating minimization (S.AltMin) to optimize (2). S.AltMin updates 𝐏 50 times in each iteration and its run-time is 50 times that of the proposed algorithm. [17] also proposes AltMin with multiple initializations for P(0) with a similarly high run-time. The results show that our algorithm (AltMin with P(0) = 𝐈 initialization) outperforms both S.AltMin and AltMin with multiple initializations. Check what AltMin refers to. Table 3. For the linked linear regression problem (Section 2.3), we report results on the three datasets from [2], available at [46][47][48], and the Ames [49] housing datasets. For all datasets, the columns of the response variables 𝐘 and the design matrix 𝐁 are centered (zeromean). 𝐁 is also replaced by its top 𝑑 principal components. For the 'scs' dataset, one of the 7 response variables is excluded to improve the model fit. For 'ftp', 'ames', 'scm', 'ames' and 'scs', the values of a feature (column of the design matrix) are rounded and data points with the same feature value are assigned to the same block and permuted.
The 'air-qty' dataset contains time-stamped readings with year, month, day, and hour information and we follow the experimental setup designed in the 'Case Study' in Section 4 of [2]: in row 4 (5) of Table 3, readings with the same month and day (day and hour) are assigned to the same block and permuted. The regression model for 'air-qty' is also as defined in [2], and we also use a moving-median filter with window size 32 to remove outliers. The proposed AltMin algorithm outperforms the competing baselines. 'AltMin-𝑘', i.e., AltMin initialized to P(0) = 𝐈, is also competitive, possibly because permuting correlated rows does not greatly corrupt the design matrix 𝐁. Results for 'Spectral' and 'Biconvex' are omitted because the methods were not competitive.
In this paper, we studied a fast alternating minimization algorithm for the unlabeled sensing problem under two structured permutation models: 𝑘-sparse and 𝑟-local. The initialization of this non-convex algorithm plays a crucial role in its performance. To address this, we proposed two initialization strategies tailored to the respective permutation models. Our primary contribution lies in the characterization of the initialization error, providing theoretical insights into its impact on algorithm performance. Additionally, we show the competitive performance of the algorithm on both synthetic and real datasets. While the current work focuses on analyzing the initialization, one direction for future research is to study the rate of convergence and establish conditions under which the algorithm converges provably to the unknown parameters of the unlabeled sensing problem.
Ahmed Ali Abbasi: Writing -review & editing, Writing -original draft, Methodology, Conceptualization. Shuchin Aeron: Writing -review & editing, Writing -original draft, Methodology. Abiy Tasissa: Writing -review & editing, Writing -original draft, Methodology. can be expressed as the difference of two Chi-square random variables. Using a tail inequality for Chi-square distributed random variables (see Lemma 4.3), Lemma B.1 provides a bound for this quantity. With this result established, the main crux of the proof of Theorem 4.5 is to relate ∑
Proof. Let 𝐄 denote the error term such that
where X(1) = 𝐁 † Ŷ(0) , Ŷ(0) = 𝐘 = 𝐏 * 𝑘 𝐁𝐗 * and 𝐄 ⟂ (𝐁) is orthogonal to the range space of 𝐁 since X(1) = argmin 𝐗 ‖ Ŷ(0) -𝐁𝐗‖ 𝐹 . Combining 𝐘 * = 𝐁𝐗 * and (39), we obtain:
The inequality in (40) can equivalently be rewritten as: 41) and applying the union bound to ∪ 𝑗=𝑚 𝑗=1 𝑗 using (47),
For 𝑡 ≥ 2 log 𝑚,
By Theorem 1.1 in [42], for a ''tall'' Gaussian matrix 𝐁 ∈ R 𝑛×𝑑 , where 𝐁 is considered tall if the aspect ratio 𝜆 = 𝑑∕𝑛 satisfies 𝜆 < 𝜆 0 for some sufficiently small constant 𝜆 0 > 0, we have Pr [𝜎 min (𝐁) ≤ 𝑐 √ 𝑛] ≤ 𝑒 -𝑐 𝑛 (see Subsection 1.2 of [42]). We consider the union bound of the failure probabilities 𝑒 -𝑐 𝑛 and 7𝑒 -𝑡∕2 . A minor calculation yields that 𝑒 -𝑐 𝑛 + 7𝑒 -𝑡∕2 ≤ 8𝑒
-𝑡∕2 for 𝑡 ≤ 𝑐 𝑛. Consequently, for 2 log(𝑚) ≤ 𝑡 ≤ 𝑐 𝑛, we obtain Pr [ ‖𝐗 * -X(1) ‖ 2 𝐹 ‖𝐗 * ‖ 2 𝐹 -2 𝑐 2 𝑘 𝑛 ≥ 𝐶 √ 𝑘 𝑛 𝑡 ] ≤ 8𝑒 -𝑡∕2 . (44) A.A. Abbasi et al. We now make a change of variables 𝑡 ′ = 𝐶 √ 𝑘 𝑛 𝑡. For 2𝐶 √ 𝑘 log(𝑚) 𝑛 ≤ 𝑡 ′ ≤ 𝐶 𝑐 √ 𝑘, we have Pr [ ‖𝐗 * -X(1) ‖ 2 𝐹 ‖𝐗 * ‖ 2 𝐹 -2 𝑐 2 𝑘 𝑛 ≥ 𝑡 ′ ] ≤ 8 exp ( -𝑡 ′ 𝑛 2𝐶 √ 𝑘 ) . (45) We make a further change of variable with 𝑡 ′ = 𝜖 𝑘 𝑛 . For 2𝐶 log(𝑚) √ 𝑘 ≤ 𝜖 ≤ 𝐶 𝑐 𝑛 √ 𝑘 , we obtain Pr [ ‖𝐗 * -X(1) ‖ 2 𝐹 ‖𝐗 * ‖ 2 𝐹 ≥ ( 2 𝑐 2 + 𝜖 ) 𝑘 𝑛 ] ≤ 8 exp ( -𝜖 √ 𝑘 2𝐶 ) . (46) (46) provides a useful (< 1) relative error bound when 𝑘 grows slowly with 𝑛, for example 𝑘 = 𝑛 𝛽 , 𝛽 < 1. □ Lemma B.1. Let 𝐏 * 𝑘 be the fixed unknown 𝑘-sparse permutation matrix with ⟨𝐈, 𝐏 * 𝑘 ⟩ = 𝑛 -𝑘. Consider the exact unlabeled sensing problem. Let 𝐱 * be the fixed unknown vector, 𝐲 * = 𝐁𝐱 * , ŷ(0) = 𝐏 * 𝑘 𝐁𝐱 * . Assuming Gaussian 𝐁, for 𝑘 ∈ [𝑛] and 𝑡 ≥ 0, we have Pr [ ‖𝐲 * -ŷ(0) ‖ 2 2 ≥ ( 2𝑘 + 4(1 + √ 3) √ 𝑘𝑡 + 10𝑡 ) ‖𝐱 * ‖ 2 2 ] ≤ 7𝑒 -𝑡 . (47) Proof. Under the 𝑘-sparse assumption on 𝐏 * , the known vector 𝐲 has 𝑘 shuffled entries. Assuming, without loss of generality, that the first 𝑘 rows of 𝐁 * are shuffled, ‖𝐲 * -ŷ(0) ‖ 2 2 ‖𝐱 * ‖ 2 2 = 1 ‖𝐱 * ‖ 2 2 𝑖=𝑘 ∑ 𝑖=1 (𝐛 ⊤ 𝑖 𝐱 * -𝐛 ⊤ 𝐏(𝑖) 𝐱 * ) 2 = 2 𝑖=𝑘 ∑ 𝑖=1 (𝐛 ⊤ 𝑖 𝐱 * ) 2 ‖𝐱 * ‖ 2 2 ⏟⏞⏞⏞⏞⏟⏞⏞⏞⏞⏟ ≜𝑇 1 -2 𝑖=𝑘 ∑ 𝑖=1 𝐛 ⊤ 𝑖 𝐱 * 𝐛 ⊤ 𝐏(𝑖) 𝐱 * ‖𝐱 * ‖ 2 2 ⏟⏞⏞⏞⏞⏞⏞⏞⏞⏟⏞⏞⏞⏞⏞⏞⏞⏞⏟ ≜𝑇 2
𝑇 1 , defined in (48), is the sum of 𝑘 independent Chi-square random variables and is bounded, using (11), as follows:
The product random variables in 𝑇 2 are distributed as the difference of two independent 𝜒
2 random variables. 𝐛 ⊤ 𝑖 𝐱 * 𝐛 ⊤ 𝐏(𝑖) 𝐱 * ‖𝐱 * ‖ 2 2 ∼ 1 2 𝑍 1 𝑖 -1 2 𝑍 2 𝑖 ∀ 𝑖 ∈ 𝑛 -𝑘 + 1, … , 𝑛. (50) The random variables (rv) in (50) are not mutually independent, but each rv depends on, at most, two other rvs. To see this, let permutation 𝐏 such that 𝐏(𝑖) ↦ 𝑗, then 𝐛 ⊤ 𝑖 𝐱 * 𝐛 ⊤ 𝑗 𝐱 * is not independent of 𝐛 ⊤ 𝑗 𝐱 * 𝐛 ⊤ 𝐏(𝑗) 𝐱 * , 𝐛 ⊤ 𝐏 ⊤ (𝑖) 𝐱 * 𝐛 ⊤ 𝑖 𝐱 * . (51) The 𝑘 rvs in (50) can therefore be partitioned into three sets 𝑃 , 𝑄, 𝑅 such that the rvs within each set are independent. Let 𝑘 1 be the number of rvs in set 𝑃 . The sum 𝑇 𝑃 , where 𝑇 𝑃 ≜ 1 ‖𝐱 * ‖ 2 2 ∑ 𝑖∈𝑃 𝐛 ⊤ 𝑖 𝐱 * 𝐛 ⊤ 𝐏(𝑖) 𝐱 * = 1 2 ∑ 𝑖∈𝑃 𝑍 1 𝑖 -1 2 ∑ 𝑖∈𝑃 𝑍 2 𝑖 , (52) is upper bounded in probability as Pr [𝑇 𝑃 ≤ -2 √ 𝑘 1 𝑡 -𝑡] ≤ 2 exp(-𝑡). (53) (53) follows from applying the union bound to probabilities 𝑝 1 , 𝑝 2 which are given by 𝑝 1 = Pr [∑ 𝑖∈𝑃 𝑍 1 𝑖 ≤ 𝑘 1 -2 √ 𝑘 1 𝑡 ] ≤ 𝑒 -𝑡 , (54) 𝑝 2 = Pr [∑ 𝑖∈𝑃 𝑍 2 𝑖 ≥ 𝑘 1 + 2 √ 𝑘 1 𝑡 + 2𝑡 ] ≤ 𝑒 -𝑡 , (55) and bounding 𝑝 1 , 𝑝 2 using the tail inequalities (10) and (11), respectively. Defining 𝑇 𝑄 , 𝑇 𝑅 similarly to (52), with cardinalities 𝑘 2 , 𝑘 3 and applying the union bound as in (54), (55) gives Pr [ 𝑇 2 ≤ -2 √ 3𝑘𝑡 -3𝑡 ] ≤ 6𝑒 -𝑡 . (57) Applying the union bound to (49), (57) gives the result in (47) with 𝐶 = 4(1 + √ 3). □