Hierarchical Multi-Stage Recovery Framework for Kronecker Compressed Sensing
Authors: Yanbin He
Geethu Joseph
Venue: ICLR 2026 Poster (2026)
Keywords: Compressed sensing
Kronecker product
Restricted isometry property
Hierarchical sparsity
Tensor operation
Primary area: probabilistic methods (Bayesian methods, variational inference, sampling, UQ, etc.)
TL;DR: We study the Kronecker compressed sensing where we demonstrate a hierarchical view of the Kronecker compressed sensing and develop a versatile algorithmic and theoretical framework.

Abstract: In this paper, we study the Kronecker compressed sensing problem, which focuses on recovering sparse vectors using linear measurements obtained using the Kronecker product of two or more matrices. We first introduce the hierarchical view of the Kronecker compressed sensing, showing that the Kronecker product measurement matrix probes the sparse vector from different levels, following a block-wise and hierarchical structure. Leveraging this insight, we develop a versatile multi-stage sparse recovery algorithmic framework and tailor it to three different sparsity models: standard, hierarchical, and Kronecker-supported. We further analyze the restricted isometry property of Kronecker product matrices under different sparsity models, and provide theoretical recovery guarantees for our multi-stage algorithm. Simulations demonstrate that our method achieves comparable recovery performance to other state-of-the-art techniques while substantially reducing run time owing to the hierarchical, multi-stage recovery process.

Published as a conference paper at ICLR 2026
HIERARCHICAL MULTI-STAGE RECOVERY FRAME- WORK FOR KRONECKER COMPRESSED SENSING
Yanbin He, Geethu Joseph Signal Processing Systems Group, Delft University of Technology ( he-he, g, j.jospphh tudel ft.nl
ABSTRACT
In this paper, we study the Kronecker compressed sensing problem, which focuses on recovering sparse vectors using linear measurements obtained using the Kro- necker product of two or more matrices. We first introduce the hierarchical view of the Kronecker compressed sensing, showing that the Kronecker product measure- ment matrix probes the sparse vector from different levels, following block-wise and hierarchical structure. Leveraging this insight, we develop versatile multi- stage sparse recovery algorithmic framework and tailor it to three different spar- sity models: standard, hierarchical, and Kronecker-supported. We further analyze the restricted isometry property of Kronecker product matrices under different sparsity models, and provide theoretical recovery guarantees for our multi-stage algorithm. Simulations demonstrate that our method achieves comparable recov- ery performance to other state-of-the-art techniques while substantially reducing runtime owing to the hierarchical, multi-stage recovery process.
INTRODUCTION
Kronecker compressed sensing (KCS) is a measurement framework that employs the Kronecker product of multiple factor matrices as a measurement matrix, capturing multidimensional signal structure while reducing measurement complexity. It appears in many acquisition systems, such as sensor arrays in communication systems (He & Joseph (025a) or separable filters imaging (Friedland et al. 2014). We focus on the general KCS problem with canonical form,
y = Hæ + n (H₁ & H1-1 @ ... Hi) x n HH) x +n.
(I)
Here, E & is the unknown sparse vector and y E R is the noisy measurements via known measurement matrix = , where each factor matrix H; E M.N. has hall row rank.
A key challenge in solving Equation is the high dimensionality of the multidimensional signal ac. grows rapidly with both the number and size of factor matrices H,, c.g., O(N') if = O(N). Another challenge is exploiting sparsity patterns as prior knowledge. Beyond simple sparsity, the nonzero elements in æ often exhibit more complex but regulated patterns. We consider three preva- lent models. The first model is the standard sparsity, where the nonzero entries be positioned arbitrarily. This model is ubiquitous and has been applied to various fields, such as image processing (Duarte Baraniuk 2010 Li & Bernal 2017 Zhao al. 2019), system identification (Sun al.
c.g.,
2022 Yuan et al. 2019). regression (Ament & Gomes 2021), and communications (Berger et al,
2010 Xiao et al. 2024). The second model, hierarchical sparsity, considers a vector partitioned
into blocks at multiple levels with sparsity structured across these levels. For example, in massive machine-type communication (Wunder al., [017]] Roth et al. 2018 2020). only subset of devices are active (device-level sparsity), and each active device sends a sparse signal, forming two-level hierarchical structured sparsity pattern on x. The third model, Kronecker-supported sparsity (or
block tensor sparsity) (He & Joseph| 2025a 2023 Caiafa & Cichocki 2013 Zhao al.) 2019 Boyer& & Haardi 2016). assumes the support x is the Kronecker product of multiple binary sup- port vectors. This pattern arises in radar imaging and wireless communications, where signals
separable across dimensions (He & Joseph| 2023) Xu et al., 2022 He & Joseph (202d)). Motivated by varied sparsity patterns, we focus efficient methods for KCS with structured sparsity.

---

Published as a conference paper at ICLR 2026
This paper introduces novel hierarchical view on KCS, showing how its dimension-win measuring structure can be used to design and analyze efficient recovery methods to exploit structured sparsity effectively. Our main contributions
follows:
We establish that when measuring via Kronecker product matrices, each factor matrix in the Kronecker product captures the vector distinct hierarchical level. It provides unified perspective for handling different sparsity models within a single framework.
Hierarchical View:
Unified Algorithm: We design a multi-stage sparse recovery algorithm using the hierarchical view. By leveraging the Kronecker structure of H through tensor operation and investigating the un- derlying structure, our method achieves a significant complexity reduction, c.g., reducing from O ((MN))) (He & Joseph, 2025a) to ((MN)) regarding stooeccer-ssppooter sparse vector covery, and accommodates the mentioned sparsity patterns within a single, flexible framework.
Theoretical Guarantees: We establish a unified restricted isometry property (RIP) analysis for KCS covering the standard, hierarchical, and stonecker-supporter sparsity. It proves that sparsity at each hierarchical level, rather than total sparsity, drives the recovery. Our result improves the RIP-based bound for KCS with standard sparsity and provides cohesive understanding of structured sparsity. We also provide a RIP-based recovery guarantee for our unified algorithm.
Related works: The Kronecker product measurement matrix is introduced for compressed imaging in in/Rivenson & Stern (009). KCS is formalized in Duarte & Baraniuk (201 Ta) tailored to hyperspec- tral imaging, with RIP analysis for KCS with standard sparsity (Duarte & Baraniuk] 201 1a b). It bounds the restricted isometry constant (RIC) of the Kronecker product using the RIC of fac- tor matrices H However, the recovery algorithm fails to leverage the Kronecker structure in H. To leverage this structure, Kronecker orthogonal matching pursuit (KroOMP) (Caiafa & Cichocki] ) adopts tensor operations. Nonetheless, it still incurs a high complexity of O(N"). ). and lacks theoretical analysis. FFieelladd et al. (214] 2015) presents two algorithms: one uses tensor unfold- ing for sequential recovery in dimension, and the other uses approximate Tucker decomposition to recover along each dimension for compressible image and video representation and recovery. Still, both approaches are limited to standard sparsity. Li & Bernal (2017) decomposes the unfolding- based approach into multiple independent subproblems for hyperspectral imaging. Yet, it fails to exploit joint sparsity patterns and is not immediately extendable to other sparsity patterns.
KCS with structured sparse recovery is also investigated in the literature. For hierarchically sparse vectors, Roth et al. (2020) discusses the hicrarchical hard thresholding pursuit (HiHTP), adapting classic hard thresholding pursuit (HTP) with tailored RIP and coherence analysis for channel esti- mation for massive multiple-input multiple-outpil systems (Wunder et al., 2019). However, it fails incorporate the Kronecker structure in H. leading to higher computational costs. For Kronecker- supported sparsity, both greedy and Bayesian methods have been explored. An orthogonal matching pursuit (OMP)-basee algorithm offers reduced complexity (Caiafa & Cichocki 2012 Caiafa & Ci- chocki, 2013) but performs poorly in noisy settings (He & Joseph|| (205aa). Bayesian algorithms, designed for applications such as hyperspectral image processing (Zhao al., 2019) and wireless communication (He & Joseph| 2025a: Chang & Su 2021 Xu et al.) (022). use structured prior distribution. They suffer from poor generalization and high complexity (He & Joseph 20256). Be- sides, both OMP-based and Bayesian algorithms lack theoretical guarantees. Recently, He & Joseph) (2025c) provides an algorithm and RIP analysis for KCS for the I = 2 case. However, the analysis decoupled from the algorithm. Also, it relies heavily on specific matrix properties, making the generalization to higher orders (I > 2) nontrivial.
To summarize, existing approaches reveal several literature gaps. First, KCS methods mostly ignore the structures of H. relying on generic solvers, while our method is specifically designed to leverage the Kronecker structure through tensor operations. Second, current methods are largely tailored to a single sparsity pattern and cannot be generalized, whereas our work provides a unified framework for multiple patterns. Third, many methods suffer from high computational complexity, while our approach is efficient and low-complexity. Besides, no prior work offers a unified RIP analysis of Fromecker-structurr matrices across various sparsity patterns, nor a recovery framework for different sparsity patterns with Rased guarantees, which are our central theoretical contributions.
Notation and tensor preliminaries: We use [4] denote the set {1, .....,,) I} for any scalar / and IN to denote the N x N identity matrix. The symbols & and x, denote Kronecker and jth mode
product, respectively. The jth mode unfolding Tú of tensor E R N1 N2 x N2 is [To]] Inj.k
2

---

Published as a conference paper at ICLR 2026
for j E [4] with = =+=++ Ei= (Ih=i Nop (ni - 1), with ny E [Ng]] Also, (oyle is (ny. k)th matrix entry, and [T] Ma, is the (M1,... ar)tth tensor entry. The ith mode product of D; E RNxxM, with T is = X, D; E IR (M.xx-xMM. X. N x.. XMM, The
ith mode unfolding of M is ( ()) (Kolda & Bader 2009).
2 HIERARCHICAL VIEW OF THE ARONECKER-STRUCTURE MEASURING
Our hierarchical view builds on the Kronecker structure in Equation interpreting the measure- ment matrix as probing the signal's sparsity across multiple block-wise and hierarchical levels. To illustrate this, we first introduce the hierarchical block partition of a sparse vector æ E R
Hierarchical partition: We first partition æ in Equation Tinto N, equal-length blocks, denoting the
Ith level blocks ((F())) my== I N. Each (((r) is further partitioned into N1-1 blocks, denoted as (I 1)th level blocks 1,n)) b 1= RIT N; We continue until reach blocks of length Ni at the second level. The first-level blocks are the individual entries of æ.
Ni-1
En to denote a block in the jth level with length IP N₁ and encapsulation nj == (M...... #T-1, 12). An encapsulation nj == (Mj, -1 can be viewed as coordi- nate for blocks in this hierarchical block structure. Also, set [[a.n]] contains all N, child blocks that share the same parent block at the level j + as that of zny- We illustrate a hierarchical partition for
For brevity,
a E IR in Figure
where [ [1.3]] = [2(2,))] = (1.3) (22,3) J as they share the parent (3)-
(a)
(1)
(2)
13
(4)
(5)
3
3.0
2 112)
2
((
(23
(2.4)
51
(2.55)
2 (1.11)) 322.2.1
41222 (3.113))
(02.2.1
3 (12.5)
(6)
X
Figure 1: (a) Hierarchical partition for 1 with N3 5, N₂ N1 4, and N 40. (b) Reordered tensor X. (c) Mode unfolding ((2) and the relation between the n2 1 row within the nyth column block and the (T 1)th level child block Fn,-1 with 74-1 = (2.1) and 3.
Hierarchical view: We first focus on the noiseless version of Equation Tefoormlatedd using tensors,
T = Y X H1 X HI, where the first mode unfolding satisfy T on the /th mode leads to
and vee vec(T()) vec(Y(n)) = y. Unfolding
To HIX (( H,U, E MIXITTI M.
R NI M; and X()) E MIXIlI N; whose nyth row the Ith level block 2n. with n, (ny). Therefore, matrix H, acts on U, and zero row in U, indicates that the corresponding /th level block is entirely zero. Hence, matrix H, captures the sparsity pattern of the /th-level blocks.
Here, =1= Xin (SL,-
For the (T 1) level, we fold U, into new tensor T, whose /th mode unfolding Túp U,, as
T x1 Hi XI-1 Hi-1 X IN..
Unfolding along its (1- - 1)th mode gives
Tu-1) H1-1X(1-1) (IN, @ (8L H M 1xN, 1xN, M;
Here, X(1-1) E -I-NI-INNNNN-- has N, column blocks, with nyth block corresponding to oa with
nr (n.). Within the nyth column block, the nt-1th row is the (T - 1)th level child block 2n
3

---

Published as a conference paper at ICLR 2026
with R4-1 (ng-1, nz). as illustrated in Figure C. The Kronecker product IN, @ ( 1 H is a block matrix and preserves the column block structure in -J Column blocks of i-1 are associated with I-level blocks, and the rows of a column block correspond to the (I - level blocks.
Hence, the zero rows in each column block of 1- indicate that the corresponding (7 - 1 1th level blocks entirely zero. Therefore, H1-1 captures the sparsity pattern at the - 1 1ehe-eee blocks.
For general jth level, we define T = X 11 Hi X2 H2 ... X, H2 Xj+1 NN ... X IN,, and
To = H) (TTRIIM
& (o))) T TRRM, I M4,
is its jth unfolding. Similar to the column block structure (I 1)th level, we have the following. Lemma 1. Consider a sparse tensor X reordered from a sparse vector x such that vec( (()) a. For the jth mode unfolding of X, i.e., X()) and with full row rank H, 's the matrix
=x X (TTILLA, @ (0))- E N, IT-ay,
can be divided into [| y N₂ column blocks. Each block is indexed by an encapsulation nj-1 with nj = (nj+1,'... ,T) forns E for k j+1,..... I. The number of nonzero rows in column block indexed by 1j+1 equals the number of nonzero blocks in [xn.]] with nj = (mg. Ri+1,"" ... nr)-
Lemma I Timplies that matrix H, actually captures the sparsity at the jth level blocks, which refer to as the hierarchical view of KCS. The above perspective can also be interpreted directly from Equation [ The Kronecker product matrix H has recursive column-block structure: each block of columns is obtained by taking the Kronecker product of a column of H, with H which itself has column block structure. This recursive structure aligns with the hierarchical partition block of æ. Hence, in this hierarchical framework, factor matrices {H:].- operate at different levels: for any p, q with p > q. Hg first measures each ęth level block of æ. the resulting measurements of all blocks are then processed by Hp. which captures sparsity at a higher level.
3 MULTI-STAGE SPARSE RECOVERY ALGORITHM
We aim to recover x: in Equation 1 from noisy measurement y. given (H, Guided by the hier- archical view in Section 2 we next present recovery framework that handles each H, sequentially. We formally define the tollowing three considered sparsity models.
Sparsity (Standard sparsity). A vector E R.N is sparse if x: contains at most s nonzeros.
Sparsity 2 (Hierarchical sparsity). A vector E N is s hierarchically sparse with 8 (81, $1-1,"'. .81) if it has a hierarchical partition defined by (NzY_= and at each level j E [7]. every set [=my]] contains at most 8; nonzero blocks.
Sparsity (Kronecker-suppooted) sparsity). A vector E N is Kronecker supported sparse if its support is the Kronecker product of 8, sparse support vectors bj E (0,1)% for j E (T).
We note that the Kronecker-supported sparsity is special case of hicrarchical sparsity, where at each level j E [7) the sj nonzero blocks En, share the same support. See Appendix B for illustration.
Our framework first solves for U, = X(n H.) from unfolding along /th mode using T(n:= H,U, + N(n-
Here, U, exhibits a row sparsity pattern where a zero row U, corresponds to an all-zero /th level block nn Thus, recovering U, from Equation is multiple measurement vector (MMV) problem and solved using MMV algorithms such as simultancous OMP (SOMP), simultancous iterative hard thresholding (SIHT), simultancous HTP (SHTP), or MMV sparse Bayesian learning (MMV-SBL).
Let the estimate of U, be Ü, with error E, modeling the estimation error and residual noise, Ú, U, E,. In the second step, we treat U, as the noisy measurement and E, as noise, reorder them into tensor T and N such that T (z) Ü, and N₂2) = E,, to obtain T - X x H X1-1 H Xr
IN, + N. Unfolding T along its (I - 1) 1)th mode as
T T(-1) H, H N(--)-
(3)
4

---

Published as a conference paper at ICLR 2026
For standard and hierarchical sparsity models, the supports of different (I - 1)th level blocks of are different. By Lemma] [ zero (I- level blocks leads to the zero rows in each column block in
Ur-1, making it a concatenation of N, row sparse matrices U for b (n1) and nı [Nr]. We thus partition Equation into N₂ independent MMV problems as (Ta-i))m = H 1un + [N(r-1))m
(1-l (
and solve them (sequentially or in parallel) using MMV solvers. Concatenating estimates Ú == [[U,-1]1. [U1-122. ... [U,-,,x,] gives the final solution, where [Ū1-1lm = [U1-1]m + E However, for the Kronecker-supporte sparsity, Equation is single MMV problem because the support is common across the (7 - 1)th level blocks.
Generalizing, for jth mode unfolding step, with measurement U,+1 from the previous step,
Úن+1 = U + Ej+1 ( (TREE @ H.)) + Ex+1-
(4)
We unfold the measurement tensor formed from Uن+1 along its jth mode as
Tu) H,U, N
(5)
Lemma 1 reduces Equation to i==f N₂ independent MMV problems for standard and hierarchi-
cal sparsity. Sparsity varies across MMVs for the standard model (defined by total sparsity rather than level-wise sparsity) but remains identical in the hierarchical model. For Kronecker-supported sparsity. Equation is a single MMV due to shared block support. While mixed models with single and multiple MMVs at different levels are possible, focus on these three main cases for brevity, leading to the Multi-Stage Recovery (MSR) algorithm, summarized in Algorithm !
Algorithm 1 Multi-Stage Recovery (MSR)
Input: Measurement y, dictionaries (He) H= E RRM.N.;
1: Fold y Y according the dimensions of dictionaries (H:) and initialize T Y
2: for j = 1, 1,... do
3:
Obtain the jth mode unfolding of i.e., Tú)
4: Solve Equation 5 for U, via compressed sensing algorithm to get estimate Ü,
5: Fold U, back to T such that the jth mode unfolding of T, i.e., T(y) is Ü,
6: end for Output: Estimated sparse vector a vec(Ü1)
Complexity:
We compare the complexity of MSR variants with existing methods for each sparsity
model, assuming Equation 5] solved sequentially, and M, = O(M), N; O(N) ) for i E [4] with I< M < N. For standard sparsity, MSR with OMP matches the time complexity of KroOMP Ca- iafa & Cichocki 2013). but reduces space complexity from O(N') ) to O O(M'~IN) For hierarchical sparsity, our MSR with HTP has time complexity O and space complexity O(MI-1N) im- proving over HiHTP (Roth al. 2020) with time and space complexities of O (22 for / 2. For -ropecker-sppported sparsity, MSR with SBL lowers time complexity to O'MN') and space complexity to (NN') compared AM- and SVD-KroSBL (He & Joseph, 2025a) with both com- plexities O(M' N The improvements are due to i) the exploitation of the Kronecker structure through tensor operation, reducing the dimensionality; and di) leveraging the MMV structure from We refer to Table in Appendix F for comprehensive comparison.
Lemma I
4 UNIFIED ANALYSIS FOR STRUCTURED SPARSITY MODELS
We establish a unified RIP analysis via a generalized notion of RIP called the (s, N)-RIP condition with : (sI, $1-1 .1) and N (N,, Ny.1. N2) defined by the dimension of factor matrices in KCS. To this end, we introduce the generalized (s, N) sparsity model, tailored to the KCS problem, which reflects hierarchical view where sparsity each level affects recovery.
Sparsity 4 (Generalized sparsity). Consider KCS with H; R M;xN, A vector R is (s, sparse if for tensor X E IRN x- reordered from x using N == (Nr, N1 .Ni,, the maxi- number nonzero rows of each of the column blocks of its jth mode unfolding X (j) is 85-
5

---

Published as a conference paper at ICLR 2026
Relation to other models: We relate the above model to the standard, hierarchical, Kronecker- supported, and block sparsity models. The standard sparsity model is not special case of N) sparsity, but the set of s sparse vectors is contained in union of (s.N) sparse vectors.
Lemma 2. Let set S contains all s standard sparse vectors in Rº, and Sa contains all Ns
sparse vectors in R for given Ns. Then, S C UsEfe(()) where fx(s) = {s E/=
+ (T - 1): 1 S 84 S 8)
Hicrarchical sparsity is a special case of (s, N) sparsity when the hierarchical partition structure matches the dimensions of factor matrices in the Kronecker measurement matrix. If, additionally, all the column blocks of jth mode unfolding (i) share the same support regarding nonzero rows, then we arrive at the Kronecker-supported sparsity. Block sparsity can also viewed as (s, sparsity with I = 2 when the block boundary matches the hierarchical partition structure.
We next define the (s, -I-PIP condition for a Kronecker product matrix H.
Definition ((s, ()-PP). A Kronecker product matrix H , with H; RRM,N, satisfies (8, RMP if there exists E (0, 1) such that for all (s, sparse x E RN, satisfies [1=8]|==13 S Hæ 1 & |)||2 [2 The smallest feasible feasible 8, denoted as 8(a,N) ()), is the (s. of H.
Under our models, (s,, N)- is defined over the unions of subspaces, thus can be used to guarantee the success of recovery algorithms, such as iterative hard thresholding (IHT) and HTP (Blumensath) 2011]. In general, such guarantees are established using the upper bound of the RICs. Therefore, we first derive the upper bound of (),, then discuss its implications for different sparsity models, and finally discuss the associated recovery algorithms and guarantees. Here, we denote the standard s-RIC of matrix H as 8, (H).
Theorem 1. The (s, N)-RIC of Kronecker product dictionary H = 0 H i.e, N) (H), satis- fies (() H)) ,1 + (H.))
The above result immediately applics to hicrarchical and Kronecker-supporte sparsity, both special cases of (s,N) sparsity. For Kronecker Kronecker-supporteed sparsity, a tighter bound could be expected due to its additional joint sparsity structure arising from the shared support across the nonzero block. However, improving the RIC bound by exploiting this additional joint sparsity is difficult. As noted in Li & Petropulu (013); Eldar & Mishali (2009), RIP analysis considers the worst-case perfor- mance and does not guarantee that MMV outperforms the SMV case. So, our bound shows no improvement, and deriving stronger RIP-based condition for the MMV model is open problem.
Theorem can also be tailored to standard sparsity using Lemma 2
Corollary 1. Consider the Kronecker product H = H For any s, the s-RIC of H satisfies (H) & maX+efu(s)
(H)
(H) I maxx I1,, + 5s, (H,)) - 1.
The s-RIC bound corroborates that only the sparsity level at different levels of blocks explic- itly affects the s-RIC of Kronecker-structured H. Also, known upper RIC bound is 8, (H) s
Il 6H,,))
Duarte & Baraniuk 2011ap. Our bound slightly improves this bound:
II + I (H.) - I 5, (,H))) - 1,
max (s)
because 8 is non-decreasing function of (Foucart & Rauhut 2013) and sf s for all i E [4] and the equality cannot be achieved simultaneously.
Maximum sparsity level: Corollary | Iiniicates that recovering s standard sparse vectors via KCS with M, < N₁ is only guaranteed when < min, N₁, as it is a worst-case analysis. When s = arg min, N,, worst-case scenario is 8j Then, 8., || H,, IN, ll 1, making H, is non-injective map, and recovery is impossible. This also indicates that it is only possible to recover block-sparse vectors with block length smaller than min, N,. However, recovery is still possible for min; N, in structured sparsity settings.
min; No with
N, and s; 1 for all j.
Measurement bounds for classical methods: We discuss the implications of Theorem on mea- surement bounds for recovering (s, N)-sparse vectors using classical iterative algorithms, namely

---

Published as a conference paper at ICLR 2026
IHT and HTP. For both algorithms, at iteration k. the support is updated via thresholding operator Lg as Th+1 = Ls (x* + H (y-Hx")) The thresholding operator depends on the sparsity model. For standard s sparse, Ls returns the support of the largest entries of in amplitude (Fou-) cart & Rauhut 2013). For & hierarchically sparse, selects the top #1 entries within each first-level block, then recursively picks top $2,.....8 8y blocks higher levels based on the E2 norm, in Roth et al. (2020). However, finding the thresholding operator Ls for Kronecker-supporct tors is NP-hard and available in the literature. For example, when / 2, it reduces to selecting rows and columns whose intersection maximizes the squared sum, equivalent to the NP-hard max- imum weight biclique problem. A practical alternative level by & norm, then recursively sum norms across matching indices at each lower level and select the top 81-1,... 81 blocks; this the approach we use in simulations for comparison. Then, IHT applies a simple projection while HTP solves a least-squares problem on the support,
sparse vec-
first select the top sy blocks at the /th
xt (a* H yy Ha*)) T
(IHT)
x k+1
arg min |yy Hall, supp(æ) E Tk+1 R
(HTP)
where operator (+) only preserves the entries within the set Th+1 and sets the others to zero.
We next discuss the implications for measurement bounds. It is known that for IHT and HTP to recover vector from union of subspaces, tailoring the thresholding operator Ls to the union and having an RIC below 1/v33 over that union is sufficient to guarantee convergence to the ground truth (Foucart & Rauhut 2013 Roth et al. 2020). So, our results shows that max ((,,N) < 1/v3 (for hierarchical sparsity) are sufficient for the success of IHT and HTP. However, it does not guarantee the recovery of the s Kronecker-supportedd sparse vectors as the thresholding operator is suboptimal.
1/v 3 (for s standard sparsity) and N(
To compare the measurement bound for KCS, we consider the simplest case with / and 8 (s152) for E (s). and Gaussian factor matrices H, 's. For recovering standard sparse vec- tors, our Corollary limplies that each H, satisfies the s, RIP, requiring M, = (s,, log N. (Foucart & Rauhut 2013). So, the total measurement bound scales M ( log N₁ log No) improv- ing over the existing bound M = (s log N1 log N₂) = (s log N₁ log N₂) (Duarte & Bara- niuk 1a) In comparison, standard compressed sensing with fully unstructured Gaussian matrix requires only S2) 2 log N2) measurements, which is smaller due to greater flexibility and ran- domness in measurement. However, KCS exploits the multidimensional structure to reduce the computational complexity during recovery. For the recovery of s hierarchical sparse vectors, Corol- lary suggests measurement bound ((s1 $2 log N₂ log N₂). while a fully unstructured Gaussian malrix requires only O(81 82 log N1 + 82 log N2) (Roth et al. 2020).
Measurement bounds for our MSR: We now establish recovery guarantees for MSR with IHT and HTP using the RICs of factor matrices.
Theorem 2. Consider the sparse recovery problem, y H x tnn. n. Define tensors X and N. which are reshaped from x and n, respectively, using the dimensions of H; 's. If x is an standard
sparse vector and the factor matrices for E [1] satisfy 83e (H,) < 1/v3 for Vs E fx(s), then
the estimate 2 of x using k-iteration IHT or HTP in Algorithm 7 satisfies
1 6-1
|| x0122
MAX
Σ
ΣI 1=
I
U.).
II T
[
(Ef
B2, myy
where [U:]m-, (ammm H.)
and if x is an s hierarchically sparse vector, and the
factor matrices Ho E []] satisfy Sz (H₂) < 1/v/3, then the estimate à of x using k-iteration
IHT or HTP in Algorithm 7 satisfies
Σ 12. .ny
ΣI go [.
Ti |||
and if x is an Kronecker Sropecker-supported sparse vector, there is
- x <
II T
|Ull
+ II T.
=
i=1
7

---

Published as a conference paper at ICLR 2026
where U; ((i) (4 N. @ (8/.-, Hi)) with ai < 1, and Ta are
MSIHT: a, = v3dg., (H.): To = (1-af) /3(11 + 822. (H,))
and
1 1 ai
v /211- 02a, (Hi)) + v da, (H:)
MSHTP: 04 = 288, (H))(1 2s; (H.)): T3 (1 (1,)
(1- 820, (H,)) (1 0:))
As the number of iterations k - 0o, the error bound reduces to 71 II'= | for the stan- dard and hierarchical sparsity, and Il Ti ||||| for Kronecker-supported sparsity. So, MSIHT and MSHTP approach the true value within a constant factor of measurement noise power. Although factors Ti IT T, and IT TV suggest error propagation the algorithm proceeds from j till j = 1 and scale with the problem dimension, this amplification is not observed in practice (see prooecckr-supporrte sparsity is tighter than that for the other two mod- els because it solves a single MMV problem, resulting a collective error bound, instead of a looser bounds due to the sum of each individual MMV bound. While our MSR's ' measurement bound scales the same as classical methods due to shared requirement on the R:-RIP of H,'S 's it can have larger error from propagation, potentially requiring more iterations H,'s with smaller s;-RICs. How- ever, a key advantage of MSR is that it provides recovery guarantees for the Kronecker-supported sparsity model, unlike classical IHT and HTP-based methods.
Figure 3). The bound for
NUMERICAL EVALUATIONS
For numerical results, we combine MSR with MMV-SBL (Wipf & Rao 2007). SIHT (Blanchard jet al. 2014). SHTP (Blanchard et al. 2014), and SOMP (Tropp et al.) 2006), and the resultng al- gorithms are referred to MSSBL,, MSIHT, MSHTP, and MSOMP, respectively. Our benchmark for the standard sparsity is KroOMP (Caiafa & Cichocki] 2013). Here, we omit computationally in- tensive SBL and OMP whose results are identical KroOMP. For hierarchical sparsity, our bench- mark is the state-of- of-tte--rttttttt HiHTP (Roth et al. 2020). For Kronecker-structured support sparsity, we benchmark with the sfate-of-the-ane AM- and SVD-KroSBL. (He & Joseph 2025a). Unlike the OMP/SBL-based algorithms, the STT/HTP--base algorithms need the true sparsity level s as input.
For all three models, we set M, M, N, N, and s; for E [T]. For the standard sparsity, we opt for H = 2|,H, H. with = 2. M 64, and N = 80. The entries of H; and the nonzero entries of æ are drawn independently from the standard normal distribution. We set s = 15, and the support is randomly drawn from a uniform distribution. For s hierarchically sparse vectors, we also opt for 1 = 2, M = 64, N 80, and s 15. Here, supports are generated by first selecting blocks uniformly at random, then assigning support within each block uniformly. In the Stonecker-supporte
sparsity model, we opt for = 3. M = 15, N = 18, and = 4. The measurement noise is zero mean white Gaussian noise whose variance is determined by SNR (dB) = 10 10 [|] H / ///m>//// |m]3]] of [3, 5, ... 23, 25).
Our metrics runtime and the normalized squared error NSE || a 12/ x L.3 where is the ground truth and & is the estimated vector. The results are shown in Figure 2 and Table T with the figure showing median and 25%/75% quartiles, and the table showing averages. The NSE for recovering an standard sparse vector is shown in Figure 2a Compared to KroOMP, MSOMP provides similar performance regarding NSE but needs one to three orders less runtime, as in Table I MSSBL outperforms KroOMP in all SNR cases with one or two orders less runtime. The NSE for hierarchical sparsity is shown in Figure 2b using only the HTP/IHT-bard parison in Appendix [H). Our MSHTP/MSIHT offers similar performance to HTP and HiHTP, and IHT. However, MSHTP requires two orders less runtime than HTP and one order less runtime than HiHTP: and MSIHT requires two orders less runtime than IHT. The NSE for Kronecker-suppor sparsity is shown in Figure 2c] Our MSSBL consistently achieves a comparable NSE and is two or three orders faster than AM- and SVD-KroSBL. In summary, MSR variants achieve similar or better accuracy than existing methods while drastically reducing computation time.
algorithms (full com-
Figure 3 shows how the NSE and runtime (median with 25% 75%% quartiles) of MSSBL. MSHTP, and MSIHT scale with the problem dimension, focusing on hierarchical sparsity. We choose = and SNR as 20dB and vary N (50,60,- ... 110). that the problem dimension N N
8

---

Published as a conference paper at ICLR 2026
10
10'
NSE
20
DUSHTP SHITE VISIIT
10
ARROHP
DIV KedRr AMARIIC
SNR (dB)
(dB)
(0B)
(a) Standard sparsity.
(b) Hierarchical sparsity.
(c) Kromecker-supported sparsity.
Figure 2: NSE as function of SNR.
Table 1: Average runtime in seconds. Bold: the best result.
SNR
dR
7 dB
11 dB
15 dB
dB
23 dB
Recovery of sparse vectors
MS0MPPSeg
0.4256
0.3827
0.2204
0.0568
KroOMP Caiafa Cichocki 2013)
130,5405
108.0526
76.6942
39. 39,9844
11,5774
0.7525
MSSBL-Seg
1.8191
1.0016
0.5758
0.2218
0.1417
0.1141
Recovery heerarchically sparse vectors
MSHTP-5eg
0.0379
0,0305
0,0297
0.0247
0.0186
0,0168
HiHTP Rath et al. 2020
0.6512
0.5493
0.5204
0.5444
0.4398
0.4574
HTP
2:2436
1.3256
0.8264
0,5311
MSIII Seg
0,0500
0.0510
0.0532
0.0509
0,0450
0.0434
IHT
8.2437
8,2412
8.2554
8,2917
8,2889
8.2789
Recovery Krronccko--upporicd sparse vectors
MSOMP
0.0142
0.0041
0,0040
0,0438
0.0026
0,0015
MSSBL
0.0728
0.0587
0.0447
0.0279
0.0119
0.0051
WVD-KreSRI, (He Jnseph 2025a
1233
5.4067
4.0681
AM KroSBL (el & Joseph 2025a)
55.9532
63.4676
75,9727
74,5840
51.7089
54.1331
125000, 216000, ................................ 1330000, where M = [(0.6N)1/17 and s [0.4N]. As expected, parallel implementation is faster than sequential. MSSBL. has the best NSE but is slower than MSIHT and MSHTP. The MSIHT is worse than MSHTP due to IHT's slow convergence (Foucart & Rauhut 2013). Overall, our MSR efficiently handles large dimensional KCS problems.
Otsiin
Damirr in Aseine QuSITIN
NEGHTF
10a
NSE
10
10
10
24
2.
2
2
Probiom Dimension
Probiom Dimension
Figure 3: NSE and runtime of MSR as functions of problem dimension N.
Application wideband massive multiple-input multiple-output (MIMO) channel estimation: Massive MIMO has been a key enabler for the fifth generation communication. For data trans- mission, an important task is estimate the channel by processing the received pilot signals sent from user. We focus on the orthogonal frequency-divise massive MIMO channel estimation, where we consider base station with a half-wavelengthh spac- ing uniform linear array equipped with Na elements serving single antenna user. Due to the environment reflection, we consider L impinging angles, each containing up to Kr delays. The maximum delay is aT, with < 1 where T, is the OFDM symbol duration. The number of sub- carriers of the OFDM symbol N₂- The channel matrix C is the superposition of impinging waves KK characterized by delays and angles as C = Σl_ Σ =i Pt. d(Ti,k, (a)) (Haghighatshoar & Caire 2017; [Cenn& & Yang 2016), where Pl, C the complex gain of the path correspond- ing to the kyth delay of the Ith angle, d(T,,))) == [L, - -.aam.t .... - (N, -1)m. /T. J is -2=(NN 1)6 T is
multiplexing (OFDM)-based wideband
the delay manifold vector of the delay Tikk, while a())) := [1, e -2200 the steering vector for O [0, representing the equivalent ith impinging angle Wunder al.

---

Published as a conference paper at ICLR 2026
Due to the significant path loss, the received signal is transmitted through limited num- ber of paths, making the channel intrinsically sparse over two sparsifying bases. The first is ob- tained by sampling the delay range [o, Tg] with Na samples as (NTa/Na)ne JN--1 leading to a de- lay basis Ha == (d()), d(Ta/Na), ... d((Na
2019).
1)],/Na]] E CNNxNa with Na = [aNNs]- The second is obtained by sampling the angular domain as {n/NVamn=0 M yielding angle basis H. = (a(0), ((a)N))... a(1 1/Na)] E CN.xN, Then the channel can be represented as C H. where X E C.NN. is the sparse representation with up to LK₁ nonzeros.
To reduce the overhead, one may use a subset OFDM subearriers and array elements for channel estimation. Denote the pilot as p CMa, where Ma < N, is the number of subcarriers in subset. Let Sa E {0, 1] Ma M,, be the sampling matrix for subcarriers and S. E {0. 1) M, Na be the sampling matrix for array where only Ma out of N, elements are chosen. We write the received signal as Y = diag(p) (Baa(p)S.CC S. + N E C XMM Mr with N being noise (Wunder et al., 2019). Plugging in the sparse channel representation and vectorizing both sides of the equation, we have
y = (() @ (diag(p)S{.t) Ha)) x n,
vec( ()(N), and (-)* is the conjugate. Denoting S, = H₂ E C.M.N,, diag(p) = H, E MaxxNa and = 2, the channel estima- tion problem is Kronecker compressed sensing problem with æ being s = (L, KL) hierarchically sparse. For simulation, we consider N₂ = 512. L {5, 10, 15, 20, 25, 35, 50, 75, 100), K2 = 3, 0.5 for the maximum delay. We fix M2
where y YY x vec(
[0.3NN] and Ma = [o.1N.], making H₂ E C154x612 and H1 E 1003x612 Both angles {0:} and delays {T are generated independently and uniformly over the sampling grid, while path gains {PI are drawn from standard normal distribution (Wunder et al.) (019). The measurement noise is zero mean white Gaussian noise whose variance is determined by SNR (dB) 10 10g E [II[((S,H)) @ (diag(p) Sa Ha)) )a: / /////a/]/// [|2]] of 20dB. We evaluate NSE := ||C H. H H P / and runtime and compare MSHTP and MSOMP to HilHT/HiHTP in WWnderr et al. 2019
N₂ 1024 OFDM subcarriers, and
Results are obtained
by two hundred independent trials.
Channel Estimation NSE
102
Hangles L
HITTF 0.0738
HIGHT
MSHIPP
MSOMP
4.1465
0.0195
0.0080
0.1303
4.1792
0.0404
0.0135
0.1444
4.1287
0.0550
0.0189
10
DMITT
20
0.1957
4.1337
0.0769
0.0447
TiilIT
0,2107
4.1670
0.0920
0.0303
MSHHTP X MSOMP
0.3117
4.1311
0.1243
0.0426
50
0.5073
4.1571
0,1816
0.0641
10-
5 10152025
35
50
100
0.8801
4.2211
0.2706
0.1888
The number of angles
1.3406
4,2272
0,2386
0.2139
Figure 4: NSE of different schemes.
Table 2: Average runtime in seconds.
We present NSE of channel estimation and average runtime in Figure and Table 2 respectively. We observe that MSHTP and MSOMP provide better performance than HiHTP and HilHT in most cases, with one or two orders less runtime. MSOMP's relatively higher NSE with large is because it wrongly identifies many insignificant paths (smaller |Pt,s, D. since it does not require the true sparsity level (L, Kr.) as input. However, we still observe that the significant paths are estimated accurately and efficiently, making MSOMP a practical option for the channel estimation task.
6 CONCLUSION
We investigated the Kronecker compressed sensing problem for signals with multiple sparsity struc- tures. We presented novel hierarchical view, comprehending that each factor matrix in the Kro- necker product dictionary senses the sparse signal at a different level, obeying hierarchical struc- ture. This insight led to a computationally efficient, multi-stage recovery framework that achieved performance comparable to state-of-the-art methods with one order or less runtime. On the theo- retical front, we unified the RIP analysis for Kronecker product matrices across various structured sparsity models, and also established the recovery guarantee for our multi- state--ttge recovery algorithm. This hierarchical framework opens promising avenues for designing new algorithms to accommo- date more structured patterns and provide efficient solutions many applications.
10

---

Published as a conference paper at ICLR 2026
7 REPRODUCIBILITY STATEMENT
All conditions required to reproduce the results are included in Section and Ap- pendix H Our implementation and data for reproducing figures and tables are available https://githb..cmm/aaaaaa Yonbennnne/nnnnn/n/m/mm
REFERENCES
Ahmed Alkhateeb, Omar EI Ayach, Geert Leus, and Robert W Heath. Channel estimation and hybrid precoding for millimeter wave cellular systems. IEEE Journal of Selected Topics in Signal Processing, 8(5):831-844 2014.
Sebastian E Ament and Carla P Gomes. Sparse Bayesian learning via stepwise regression. In Proceedings of the International Conference on Machine Learning, pp. 264-274. PMLR, 2021.
Christian R Berger, Zhaohui Wang. Jianzhong Huang, and Shengli Zhou. Application of compres- sive sensing to sparse channel estimation. IEEE Communications Magazine, 48(11): (114411644177, 2010.
Jeffrey D Blanchard, Michael Cermak, David Hanle, and Yirong Jing. Greedy algorithms for joint sparse recovery. IEEE Transactions on Signal Processing, 62(7)::694-1704 Jan. 2014.
Thomas Blumensath. Sampling and reconstructing signals from a union of linear subspaces. IEEE Transactions Information Theory, 57(7):4660-4677 July 2011.
Rémy Boyer and Martin Haardt. Noisy compressive sampling based on block-sparse tensors: Per- formance limits and beamforming techniques. IEEE Transactions on Signal Processing. 64(23): 6075-6088 2016.
Cesar F Caiafa and Andrzej Cichocki. Block sparse representations of tensors using Kronecker bases. In Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing, pp. 7709-2712, Mar. 2012.
Cesar Caiafa and Andrzej Cichocki. Computing sparse representations of multidimensional sig- nals using Kronecker bases. Neural Computation, 25(1):186-22 220, Jan. 2013.
Wen-Che Chang and Yu T Su. Sparse Bayesian learning based tensor dictionary learning and signal recovery with application to MIMO channel estimation. IEEE Journal of Selected Topics in Signal Processing, 15(3):847-85,, Apr. 2021.
Zhilin Chen and Chenyang Yang. Pilot decontamination in wideband massive mimo systems by exploiting channel sparsity. IEEE Transactions on Wireless Communications, 15(7):5087-5100, 2016.
Marco F Duarte and Richard G Baraniuk. Kronecker compressive sensing. IEEE Transactions on Image Processing, 21(4):9945004 Aug. 2011 la.
Marco F Duarte and Richard G Baraniuk. Kronecker product matrices for compressive sensing. Rice University, Houston, Technical Report, Mar. 2011b. URL httt:://ww..eee ecs.wwwwwwcwcc.o..oo. umass. edu, induarte/ wwwaaee//ooneeeeeeeeee.
pdf
Marco F Duarte and Richard G Baraniuk. Kronecker product matrices for compressive sensing. In Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing, pp. 3650-3653, Mar. 2010.
Yonina C Eldar and Moshe Mishali. Robust recovery signals from structured union of subspaces. Oct. 2009.
IEEE Transactions on Information Theory,
Simon Foucart. Hard thresholding pursuit: algorithm for compressive sensing. SIAM Journal on numerical analysis, 49(6))2543-25633 2011.
Simon Foucart and Holger Rauhut. A mathematical introduction to compressive sensing. Birkhäuser Basel, 2013. ISBN 081 0817649476.
11

---

Published as a conference paper at ICLR 2026
Shmuel Friedland, Qun Li, and Dan Schonfeld. Compressive sensing of sparse tensors. Transactions on Image Processing, 23(10):4438 4447, Oct. 2014.
1EEE
Shmuel Friedland, Qun Li, Dan Schonfeld, and Edgar A Bernal. Two algorithms for compressed sensing of sparse tensors. Comprexsed Sensing and its Applications: MATHEON Workshop 2013, pp. 259-281. Springer, 2015.
Saeid Haghighatshoar and Giuseppe Caire. Massive mimo pilot decontamination and channel in- terpolation via wideband sparse channel estimation. IEEE Transactions on Wireless Communica- tions, 16(12):8316883)) 63332, 2017.
Yanbin He and Geethu Joseph. Structure-aware sparse Bayesian learning- based channel estimation for intelligent reflecting surface-aided MIMO. In Proceedingy of the IEEE International Confer- ence on Acoustics, Speech, and Signal Processing, pp. 1-5, Jun. 2023.
Yanbin He and Geethu Joseph. Bayesian algorithms for Kronecker-structured sparse vector recovery with application to IRS-MIMO channel estimation. IEEE Transactions on Signal Processing. 73: 142-157, 2025a.
Yanbin He and Geethu Joseph. Efficient off off-grid Bayesian parameter estimation for Kronecker- structured signals. IEEE Transuctions on Signal Processing, 6--630, 2025b.
Yanbin He and Geethu Joseph. A hierarchical view of structured sparsity in Kronecker compressive sensing. In Proceedings of the European Signal Processing Conference, pp. 1-5, Sept. 2025c.
Yanbin He and Geethu Joseph. Kronecker-structured sparse vector recovery with application to IRS- MIMO channel estimation. In Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing, pp. 1-5, Apr. 2025d.
The Math Works Inc. Matlab version: 24.2.0.292308. (r2024b), 2024.
Tamara G Kolda and Brett W Bader. Tensor decompositions and applications. SIAM Review, 51(3): 455-500, 2009.
Bo Li and Athina P Petropulu. Structured sampling of structured signals. In Proceedings of the IEEE Global Conference Signal and Information Processing, pp. 1009-1012. IEEE, Feb. 2013.
Qun Li and Edgar A Bernal. An algorithm for parallel reconstruction jointly sparse tensors with applications to hyperspectral imaging. In Proceedings of the IEFE Conference Computer Vision and Pattern Recognition Workshops, pp. 24-31, 2017.
Yuxing Lin, Shi Jin, Michail Matthaiou, and Xiaohu You. Channel estimation and user localization for IRS-assisted MIMO-OFDM systems. IEEE Transactions on Wireless Communications, 21 (4):2320-2335, Apr. 2021.
Zhendong Mao, Xiqing Liu, and Mugen Peng. Channel estimation for intelligent reflecting surface assisted massive MIMO systems-A deep learning approach. IEEE Communications Letters, 26 (4):798-802, Jan. 2022.
Yair Rivenson and Adrian Stern. Compressed imaging with a separable sensing operator. IEEE Signal Processing Letters, 16(6):449-45 452, 2009.
Ingo Roth, Axel Flinth, Richard Kueng, Jens Eisert, and Gerhard Wunder. Hierarchical restricted isometry property for Kronecker product measurements. In Proceedings of the Annual Allerton Conference on Communication, Control, and Computing, pp. 632 638. IEEE, 2018
Ingo Roth, Martin Kliesch, Axel Flinth, Gerhard Wunder, and Jens Eisert. Reliable recovery of hier- archically sparse signals for Gaussian and Kronecker product measurements, IEEE Transactions on Signal Processing, 68:4002-4016, June 2020.
Luning Sun, Daniel Huang, Hao Sun, and Jian-Xun Wang. Bayesian spline learning for equation discovery of nonlinear dynamics with quantified uncertainty. Proceedings of Advances Neural Information Processing Systems, 35:6927-6940, 2022.
12

---

Published as a conference paper at ICLR 2026
Joel A Tropp, Anna C Gilbert, and Martin J Strauss. Algorithms for simultaneous sparse approxi- mation. Part I: Greedy pursuit. Signal Processing, 86(3):572-58)) 2006.
E. van den Berg and M. P. Friedlander. Probing the pareto frontier for basis pursuit solutions. SIAM
Journal on Scientific Computing, 31(2):90-911, 2008. doi: 10.1137/0807 14488.
E. van den Berg and M. P. Friedlander. SPGLI: solver for large-scale sparse reconstruction, December 2019.
N. Vervliet, O. Debals, L. Sorber, M. Van Barel, and L. De Lathauwer. Tensorlab 3.0, Mar. 2016. Available online.
Peilan Wang, Jun Fang, Huiping Duan, and Hongbin Li. Compressed channel estimation for intel- ligent reflecting surface-assisted millimeter wave systems. IEEE Signal Processing Letters, 27: 905-909, May 2020.
David P Wipf and Bhaskar D Rao. Sparse Bayesian learning for basis selection. IEEE Transactions on Signal Processing, 52(3):1533,116, Aug. 2004.
David Wipf and Bhaskar Rao. An empirical Bayesian strategy for solving the simultancous sparse approximation problem. IEEE Transactions Signal Processing, 55(())37044311) June 2007.
Gerhard Wunder, Ingo Roth, Rick Fritschek, and Jens Eisert. Hihtp: A custom-tailored hierarchical sparse detector for massive MTC. Proceedings of the Asilomar Conference on Signals, Systems, and Computers, pp. 1929-1934. IEEE, 2017.
Gerhard Wunder, Stelios Stefanatos, Axel Flinth, Ingo Roth, and Giuseppe Caire. Low-overhead hierarchically-sparee channel estimation for multiuser wideband massive mimo. IEEE Transac- fions Wireless Communications, 18(4):21862219 2019.
Zhenyu Xiao, Songqi Cao, Lipeng Zhu, Yanming Liu, Boyu Ning, Xiang-Gen Xia, and Rui Zhang. Channel estimation for movable antenna communication systems: A framework based on com- pressed sensing. IEEE Transactions Wireless Communications, 23(9):11814-11
2024.
Xiaowen Xu, Shun Zhang, Feifei Gao, and Jiangzhou Wang. Sparse Bayesian learning based chan- nel extrapolation for RIS assisted MIMO-OFDM, IEEE Transactions on Communications, 70(8): 5498 5513, Aug. 2022.
You You, Li Zhang, Minhua Yang, Yongming Huang, Xiaohu You, and Chuan Zhang. Structured OMP for IRS-Assisted mmWave channel estimation by exploiting angular spread. IEEE Trans- actions on Vehicular Technology, 71(4):444444443 Apr. 2022.
Ye Yuan, Xiuchuan Tang, Wei Zhou, Wei Pan, Xiuting Li, Hai-Tao Zhang, Han Ding, and Jorge Goncalves. Data driven discovery of cyber physical systems. Nature Communications, 10(1): 4894, 2019.
Rongqiang Zhao, Qiang Wang, Jun Fu, and Luquan Ren. Exploiting block block-sparsity for hyperspectral Kronecker compressive sensing: tensor-based Bayesian method. IEEE Transactions on Image Processing, 19:1654-1668 Oct. 2019.
13

---

Published as a conference paper at ICLR 2026
A PROOF OF LEMMA I
The proof proceeds in two parts: first, establishing the column block structure, and second, analyzing the sparsity of each block.
For the first step, we observe that matrix (Tn= N; (H))
block-diagonal matrix.
has I N₂ identical diagonal blocks, each equal (8)=)- To match this structure, we ++1 partition the columns of the unfolded matrix X)) into I N, column blocks. The standard column ordering in tensor unfolding places elements with higher-level indices (Mj+1,.... further apart. Consequently, we can partition (() into II N₂ column blocks, where each block corresponds
unique encapsulation 1j+1 (Mi+1....,,)) n1) as
X (5) [ (f)....... ...
Since (TTRET @ (8) H.) is block-diagonal, the multiplication with X (3) decouples and operates on each of these blocks independently.
X (41ttr N, @ (H.))-
X (().......)) HO)
U (o)=)-
This confirms that the resulting matrix is also composed of IT==y+1 N; column blocks, each indexed by Ni+1 and given by (()mn (0) H.)
For the second step, consider a column block indexed by a fixed nj. X (J)) H.)) The rows of this block are indexed by n, E [N,]- The kth row of X (1)) (0) (0)=)-1 H.) will be nonzero if and only if the keth row of X (g),n contains nonzeros due to the full rankness. Moreover, the kth row of (()). is the hierarchical block In, where the encapsulation is nj = (k, ni-1o..., ,1y). This equivalence follows because the indices of the entries the kth row of X (i) align exactly with those of the hierarchical block ÆnZ with nj (k, Mj+1........ m2). Hence, the kth row of X (j))m)-= and the hierarchical block Eny with nj (k, mj+1. ................................ n2) contain identical entries with identical order. Thus, the number of non-zero rows in X (f) (ol=)- H.) is the number of hicrarchical blocks { J (within the parent block defined by n;+1) that contain at least one non-zero element, which concludes the proof.
Illustrative Example: We provide an example in Figure for the proof of Lemma 1 In Figure 5(a), we consider the same vector x E R440 with = 3, N = 5, N2 = 2, N1 = 4, as in Figure and mark 3 nonzero entries using colored cubes. Figure 5tb) (c), and (d) illustrate the reordered tensor X, its mode unfolding X()) with j 2, and how U₂ computed, respectively.
The first step of the proof corresponds to Figure ((c). To why there is block column structure, we first investigate how the unfolding matrix X (2) is obtained. Since the unfolding tensor mode j - the row of the unfolding matrix ((2) is indexed by n2 1,2. The column index k is determined by ny and n, jointly i (n, 1) + N₁ (nz 1). according to the definition in Section When ny increments by one, k increases by 1; when ng increments by one, h increases by N1.
To arrange the columns of ((2). we fix ny and let n, runs through 1, 1,2,3, 4, and then increase n2 by one and let n runs through 2, 3, 4 again, shown in Figure ((c). This indicates that ny indexes I N₂ N₂ 5 column blocks, each containing IP:11 N, N₁ columns. Besides, the matrix (I N. @ (H) in this case reduces to the block diagonal matrix (Is Hi) in Figure d), matching the column block structure of X(2) Therefore, U2 can be divided into No column blocks given by the product of the column blocks of X(2) and Hi, where each block of U2 is also indexed by an encapsulation n3 (ng) for n3 E [Na]-
For the second step, to understand why the number of nonzero rows in a column block indexed by 1i+1 in U, equals the number of nonzero blocks in [mny]] with nj (nz, (nj, "if1,' ,,)),
14

---

Published as a conference paper at ICLR 2026
examine Figure d). Consider ny 4. In U2. it corresponds to the fourth column block given by (2). HJ, where (2).(4) is the fourth column block (()) Also, we have
a2p [.,)]]] = f (1,4 (2.4)))
of the column block (2). Only (22,4) is nonzero leading nonzero row in (4)).(,;Hii This demonstrates that the number of nonzero rows in column block indexed by nj+1 in U, equals the number of nonzero blocks in [æn,]] with nj = (n,, "" ,11).
Each element of [[)[ corresponds to
(
(1)
(3)
J(4
(()
(1.1)
(2,1)
(122)
(2.2)
1133)
22,33)
(14
(2,4
(1,6)
(2.6)
2
((1111
(222.1)
(3.122) (4.22) (3.1.3)
(3,2,4)
(4.255)
(8)
(c)
2
4
5
n2 = 1, 2
122
71
1, 2,
(d)
HI)
H
U2
H
H
X(2)
X
H
(Is H1)
X (2)
Figure 5: An illustrative example of the proof of Lemma with Hi R3x4
B ILLUSTRATIONS OF DIFFERENT SPARSITY PATTERNS
In this section, we provide examples of sparsity patterns considered in this paper in Figure using the same vector shown Figure
In Figure ((a), we present the standard sparsity with = 3. Three nonzero entries Ff3.1.1)- (4.1.33) and (3,22.4) are arbitrarily positioned. Take 233.2.4) an example. Its encapsulation (3,2,4) means (3.2244) is the third entry of the block indexed by encapsulation (2, 4), i.e., ((224)) while (2244) means is the second block of the block indexed by encapsulation ())))- Then 2 is the fourth block of vector T.
In Figure b)),, we show an example s (s3, $2, 51) (2, 2) hierarchical sparsity. For the third level blocks, the set [x(2]] contains all blocks that share the same parent block as x(2)) meaning [xxe)]] = (m(1)) (12) ((3) (()) (35) = [a] = [=]
[=())]]- Since $3 = 2, according to the definition of Sparsity 2 [x] contains at most 83 = 2 nonzero blocks,
[x
which are 222) and (i4))- an)) For the second level sparsity $2, take (22) and its child blocks (1.2) and
(2.2) an example. Since 82 = 1, it means that [=<],2)] (11112) (2.2)) = [={222)] contains at most nonzero block, which is (112)). Similarly, [.]].,)] contains at most nonzero block. For the first level sparsity 81 2, we take (11)2) and its child blocks as an example. There should be at most 81 = 2 nonzero blocks in the set of children of (1.2)), which are 222.1.2) and (3, in [22]. Since this is the first level, a block corresponds to an individual element of æ.
Figure ( illustrates the s = (sg, 82, $1) = (2, 1,2) Kronecker-supported sparsity with by = 82, [0, 1, 0, 0), b2 [1,0), 0), and b 1, 1, 0). Its support is then by 60 b2 69 bj.
15

---

Published as a conference paper at ICLR 2026
(8)
2
322)
23
((4)
25
11.))
(2.1)
L (22) (12
(131
(2,3
004
(2,4)
(1. (2.5)
0
5
(11.I) (2,2.))
(3.12)
(4,22)
27 (31,33
(3,244
(4.255)
(b)
1
22
33
24
(5)
(1.D)
T 2,11
T 1122
T (2,2)
(1.3)
(2.3)
22 114
(2,4)
(1.))
(2.5)
0
(3,,2)
2 (11,4
(2,2,))
(c)
x
(1)
2222)
33)
(4)
(()
(1.1)
(211
(2,2)
2 (1.3)
F02,3)
2 (1.4
Z(2.4)
21.5)
(2.6)
1
j
(3,1222
2 (11144
(2,20,0))
Figure 6: Examples of different sparsity patterns.
(a) Standard sparsity. (b) (2, 1, hierarchical
sparsity. (c) (2, 1,2) Kronecker-suppported sparsity.
C PROOF OF LEMMA
Let æ be s-sparse vector. We denote ky as the total number of nonzero blocks within all jth- level blocks of æ. Clearly, kı = 8 and ki+1 = 1. Then, each nonzero block in the j + 1th level can have at most k; (kj+1 1) number of nonzero jth level blocks. This occurs in the most unbalanced case, where k;+1 - blocks have only one nonzero jth level block while the remaining block has hig - (li+11 1) nonzero jth level blocks. This observation leads to the upper bound for the sparsity level, 8j I k2 (kz+1 - 1), which yields El I s -1.. 1). So, any sparse vector (sE
D PROOF OF THEOREM
For any æ, note that Equation bounds |Ha| Following the hierarchical view, note
H 11 X H. X H, 2 ( (8|=-1 H.) l
Using the RIC of H,, have
(1-
11,6,))| (1) (0L H,
H 2/2 5,,) (() (01=1-1 Hi) 1
We also note that | (n) (8) & tensor folding and unfolding. Therefore, using RIC of H 1-1, arrive at
H, | = H X (IN, @ (8 || due to the
(1 5-, (1 1) (IN, (0) H.)) 2
< H < (1 +1+ 8_, (1 + 68,_,))| ([x(()) (H))-, 1|
Repeating these steps recursively, following the analysis in the hierarchical view, obtain
II (x (1) (8 (8L=111 < H
I (1 + & (0=)111 1
Since X(v) (o! IN.) = X(1) and X ||
1 I (1 5z, 1
Hx
II (1 + 5 1
16

---

Published as a conference paper at ICLR 2026
Hence, we derive
(a.N) (H) S max
II 5,,), II 8z.)
1]
II i=1
5z,) 1,
which completes the proof.
We note that the high-level proof strategy of our Theorem [ and Roth et al. (200)) The- orem is similar in that both aim to sequentially unwrap the effect of the Kronecker product. The key difference is that we employ tensor representations and operations such as tensor unfold- ing, enabling a straightforward, flip-operator-free proof. This formulation clearly demonstrates how the sparse signal æ (or its tensor form X and its unfolding (() is measured by factor ma-
Remark:
trix Hj through linear transformation X() (T1E)) N. @ (o) The row sparsity of
( (TTE=t) N; @ (8) T
is dictated by the sparsity our hierarchical block partition as in Section The aspect of this multi-stage measurement framework missing in Roth al. (2020). Thus, Roth et al. 2020 Theorem 4) focuses solely on hierarchical sparsity while our multi-stage framework provides general perspective that defines generalized sparsity, where standard, hier- archical, and Kronecker-supported sparsity are special cases for analysis and recovery. This proof also explains why standard RIP cannot be improved beyond hierarchical sparsity, clarifies the max- imum achievable sparsity level, and shows why the corresponding bounds are fundamentally tight. further provides insight into why proofs for pronecker-supporte sparsity can strengthened, drawing analogies to standard RIP and MMV analyses.
E COMPLETE RESULTS ON THE NUMBER OF MEASUREMENTS
In this section, we present the measurement bounds for unstructured H with different sparsity pat- terns. Let H E IR A has independent and identically distributed standard Gaussian. For
e.V
M = In
where e is positive constant, sparse vectors can be recovered from the measurement of H with high probability (Foucart & Rauhut 2013). Also, if
Σ11 ==
M
In
8
hierarchical sparse vectors can be recovered from the measurement of H with high probability (Roth et al. 2020). These two results lead to the discussed measurement bounds in Section
F COMPLEXITY COMPARISON
We comprehensively analyze the complexity of our MSR algorithm to demonstrate the benefit of exploiting the Kronecker structure of И via the hierarchical view. We consider MSR combined
with MMV-SBL (Wipf & Rao 2007), SIHT (Blanchard et al. 2014), SHTP (Blanchard et al. 2014). and SOMP (Tropp al. 2006) sparse recovery algorithms,, referred MSHTP, and MSOMP, respectively. We also use Seq and P to represent the sequential and parallel implementation of Equation 5 Assume M, 's are ((M), , N,'s ss are O(N) for i E [4]. and M
MSSBL, MSIHT.
N. We compare the time and space complexities of our algorithms with those of other state-of-the- art algorithms. For the recovery of s sparse vectors, we include SBI. (Wipf Rao 2004), OMP, and KroOMP (Caiafa & Cichocki 2013) as benchmarks. For the recovery of hierarchically sparse vectors, HiHTP (Roth et al.) 2020). IHT, and HTP are used as benchmarks. We note that only the exact implementation of HiHTP for I = 2 is given in (Roth et al.) 2020 Regarding recovering Kronecker-supported sparse vectors, we consider AM- and SVD-KroSBL. (He & Joseph 2025a) for benchmarking.
For the recovery of s standard sparse vectors, our MSSBL and MSOMP substantially reduce both the time and space complexity compared to their traditional counterparts. In terms of time complexity,
17

---

Published as a conference paper at ICLR 2026
our MSSBL (O for Seq and O N) for P1) is superior than SBL (O(M2) while the time complexity of MSOMP (O(MN') for Seq and O(M'N) for P1) is also lower than OMP with O(M'N') Moreover, both MSSBL and MSOMP avoid (M " in space complexity and have O(MI-IN) for Seq and ) for P1. Compared to KroOMP with time complexity O(MN') and space complexity O(N'), MSOMP-Seq achieves the same time complexity but with much lower space complexity NN). Alternatively, we can achieve a much lower time complexity (MM* N)) by parallel implementation, the cost of slightly higher space complexity of O(MN)).
The computational gains are particularly significant in the context of structured sparsity. For both hierarchically sparse and Kronecker-suppont supported sparse vectors, classical methods like IHT and HTP exhibit a time and space complexity of O (MMNN ). Our MSIHT-Seq, MSHTP-Seq, and MSSBL- Seq have time complexity O O(MN2. ), ((MN') ). and ( respectively, and ((M N) for space complexity. Compared to HiHTP, our MSSBL-Seq has the same time complexity ( N2) while MSSBL-P1 has lower space complexity (O(MN)) compared to (MM) of HiHTP.
Similarly, for Kronecker-supportedd sparse recovery, when compared to AM-KroSBL and SVD- KroSBL, the MSSBL algorithm demonstrates lower time complexity from O(M*N*) to (M)N') and space complexity from O to O(N'). MSIHT and MSHTP exhibit the same or even lower time and space complexities than MSSBL, which is lower than AM-KroSBL and SVD- KroSBI., demonstrating the superiority of our multi-stage framework. We list all the time and space complexity of the algorithms in Table 3]
We use REM, RoMr, RHTT, RurT, and RAM to denote the number of EM, OMP, HTP, IHT, and AM iterations. These values can vary for different algorithms and experimental settings.
Table 3: Complexity of different algorithms in different sparse recovery problems.
Method
Time Complexity
Space Complexity
Recovery of sparse voctors
MSSBL-Seq
(REM N' + MN'))
O(M' N)
MSSBL.PP
O [REM(TM2M
M NJ)
O(MN')
MSOMP-Seg
0 RoMP MN + ROMI NI R Ro
O(M N)
MSOMP-P1
O RoMP N RouP M' + ROMP M 1
O(MN')
KroOMP
[RoMPMN' HoMp A + Bomp MN + RoMP)
O(N')
SBL.
O "!MM
O((MN))
OMP
O (Rooo((MM (MN)' + ROMP + Ro)
(((MN))
Recovery of a hierarchically sparse vectors
MSSBL-Seg
(RRM(MNN N' + MN'))
O(M'-1X N)
MSSBL-P1
0 (RPM N + M))
O(MN') )
MSHTP-Seg
0 RM + max; M(NN-1 1))
O(M -N)
MSIHT-Seq
O (R")
O(M" NN
HIHTP et al. (000 [1 = 2)
O RRry ((1ee)(()))))))))))) " + (MN)))
O((MN))
THT
0 Rayy (M
O((NNY)
HTP
0 (Rarr (MMY + (TE 4)) M)))
O((NN))
Recovery of Kronooker-supportd sparse vectors
MSSBL
0 (REMM(IM' MN
O(N')
MSIIT
(
O(N')
MSHTP
Rgrp RgrpPNN M/
RHTE M EJ 20))
O(N')
MSOMP
0 Rozap N'11 RoMP MN ROMP MN')
ONN")
AM KroSBL He & 2025u
2025a
[RMa(RaaIN + (MNN'I)
((MNN))
SVD KroSBL Hee & Joseph (2025a) (2025)
(=7 + (MNY))
HIT
(NN')
O N
HTP
O (MMr + (ITL,, M')
((MN)))
Table 3 compares the complexity of different versions of MSR to different traditional compressed sensing algorithms. The conclusion can be extended to a more general case. Consider N M and I > 1 with M, I e Z. In a prototypical compressed sensing problem yp- = HpEp with Hp E IR XN,, being a dense, unstructured measurement matrix, consider a general compressed sensing algorithm that has complexity OMM, to recover ap with a, b 2 1. We note that a, b 2 is a fair consideration since computing Hpmp already requires O(M, Np).. Then applying this algo- If we combine the same compressed sensing algorithm with our MSR, regardless of the special structure or the MMV property of Uj. step 4 in Algorithm 1 is simply solving NI-J M compressed sensing subproblems where each
rithm to Equation I indduess a time complexity of O(M"
18

---

Published as a conference paper at ICLR 2026
has M measurements and N unknowns, inducing per level complexity of Considering all steps from j = I to 1, the total complexity is given by
M M M M(N' M²) ==1
M N'. M
No
N' - -M)
- M
To compare to Mal Nb1, we consider the ratio
MCN M')
N' NM' M Nb (N - M) M=U--INNMI-D( M)
(1) M
N' M
(6) I (si)
S
M I 1,
(N M)
where (i) holds since the ratio a decreasing function of a and b, (ti) due to
N' M = (N M) N < IN M),
and (iit) holds since M 2 2'-1 2 for IA, M 1,I, Z, conclude that N/T 1, and thus M M N6 -M N-M < N for any N M > 1, 1, and b 2 1 with M, N, E E. Hence, in general, our MSR has lower computational complexity than traditional compressed sensing algorithms when both are applied to Equation [] We note that when / Equation reduces to the traditional compressed sensing problem and MSR has identical complexity to that of a traditional compressed sensing algorithm.
G PROOF OF THEOREM
Before the proof of Theorem
we introduce four aiding lemmas.
Lemma 3. (Foucart & Rauhua 2013 Lemma 6.16) Given vector N and index set T C [N], there is
[[(IN - H H) & 6.11v112,
if the cardinality of the union of and the support set of v is not exceeding L.
Lemma 4. (Foucart Raulut 2013 Lemma 6.20) Given vector n RR and set T [N] with cardinality not exceeding s, then
H n)rll
V1 82
Lemma's 5. For sparse matrix X with row support T with card(T) < and N E R MxN the
sequence (X*} defined by SIHT or SHTP for solving an MMV problem Y HX + N with
X O, satisfies for any k 2 0,
|x* - -X11= S a² a +
where
k
for SIHT: = /300o: T 13(1 820)
1- and
282 82
(12(1-
-22) + VI + -aat
for SHTP:
(1- a1)
Proof. The proof closely follows the technique in Foucart Rauhut (2013 Theorem 6.18) and Blanchard et al.] (2014). Here, we extend the SMV case in Foucart & Rauhut (2013) Theorem 6.18) to the MMV case. In the MMV case, the thresholding operator retains the rows of X + H Y H with the largest row lg norms, and then we have
|| (x* + H (Y Hxx lf < (x* H (Y H)x 11..
Removing the common rows from both sides, we arrive at
(x* H (Y H ) 1 < x H (Y HX*)) TRINT 11.
19

---

Published as a conference paper at ICLR 2026
IHT proceeds with X (x* + H (Y -HX (X)TE O, get
T.+.. ) Since (X TiTT+
0 and
|l (x - x k+1 + x x - x + H (Y - -HX*)) TYTa+1 l|F
-xx*)) yt s (x - X - (Y-
Applying reverse triangle inequality to the left-hand side and rearranging, we arrive at x TITEE (x-x*+I
|| (x*
+ H (Y --=*))) T l|f
-
+ +11 X X H (Y HX )) T llF
< v 2 (x* X + H (Y
lIf,
where TATk+1 = (T) T (T( (T++1 7) denoting the symmetric difference of the sets and Tht1 Therefore, we obtain the error in the kth iteration as
T(T+-1 1P xh+1_ XI? = 11 x) 7u+a 11 + (X k+1
X)
= =I(x H (Y - -x**) T)111 11: +11 +++1((ktt)) (x -X)
$ 11 $1( + H (Y - (Y-HX* - -X) -X) Tk+1 11
+ 2| (x* -X H (Y- H TATE:1
- 3|| (x+ * - + H T(Y HX*) ( TUTEII I1
Considering Y HX + N, we then have
||
1|Ir v /3l (x* X + H (Y - -HX* TUT--I ||
v3ll ((I H H X -x) H N) TUT lFr
< /3l ( H) (x* -x--Xi -x)) TUTU- |F + v3I1 (H" N)
TUT I|F
S 33 - X|1 + /3(1 + 826 |N
where the last step is the direct consequence of Lemma
3 and Lemma
To see this, we note
(I H H (x* x))
1/ || H
X X xIn)
||| X
1x X 1|?
35
We derive similar argument for vI + d2s N lp. Conclusion for HTP has been given in Blanchard et al. 2014 Theorem 3). This concludes the proof.
where [X* X)n is the nthh column of matrix X
X.
Lemma 6. For the sparse recovery problem in the stage of unfolding jth (j S 1- mode of tensor
T 1 Hi X2 H2 xx H, Xj+1 IN ... X1 IN₂ +N,
where the sparse tensor X corresponds to s standard sparse x or hierarchically sparse æ. Its jth mode unfolding leads to
Te H,U,+ = -, (MRe @ (H))- + +NG
(6)
Then the estimate of U, denoted U, and obtained through IHT HTP, satisfies
[|U, [U,]on S a [||U,) (n + +T,ll(00, - Ustalanaale,
where ÜŪlmis and (Ulmm == bm (oij-1H,) H.) denote the 1i+1th column block of Ú and Ug respectively. Here, encapsulation nj+1 == (nj+1. block. The indices (ni) in encapsulation nj+2 == (1242. - n1) have the same value as the indices (n in encapsulation nj+1- i.e., the block indexed by 1j+2 should be parent block of the block indexed by ni+1- Constants μ2 and T, depend on the iteration number and matrix H,-
ni) is the index the column
20

---

Published as a conference paper at ICLR 2026
solve Equation by separating it into IT: N, MMV problem, where each MMV problem is indexed by an encapsulation n1+1- Suppose we consider a fixed
Proof. According to Lemma
encapsulation - (1j+1r'" ,n]), and consider the MMV problem indexed by nj+1
(T)laj H + (N6)) (X(X())) (n)i (H)))- + [
According to Lemma and denoting the solution as [U, (nj with IHT or HTP iterations, we have 10ا UUy0n lls S SafllUU 11 ||r + T, 1,I[[N Ille
where My and ry relate to the RICs of matrix Hy. The only step left is to bound [| (1) ||F using [B5+innjaa where E = Ū,+1 - U,+1-
We recall that Equation 6 is obtained by unfolding the mcasurement tensor formed from the matrix Ú,+1 = UJ+1+ along its jth mode. Hence, N6) is simply reordered version of E Conse-
quently, the entries of the matrix [N()), are essentially entries of the n3+1 th row of [B)+1la/4** leading to
|| |I 11
and we arrive at the desired result.
To elaborate, we first investigate the indices of the entries of [N The entries of the ni row of matrix [N are obtained by i) fixing n2 nz (row index) and n1j11 " .?s nī (encapsulation), and ii) running R1, R2-1 from one till N1. Ny-1, respectively. Thus, the entries of matrix (1)), can be obtained by i) fixing Nj+1 " n1 = n7. i1) run- ning " - from one till N1, NN-1, respectively, and idi) running nj = 1,... .N (going over all rows). Given such knowledge, we start investigating the 2 th row of matrix [Ejtilnjaa The entries of this row are obtained by i) fixing Mj+1 = 25+1 (row index), i1) fix- ing nj42 "j+2.'..,, nr ni (fixed encapsulation), and wi) running 71..... ,, from one till N,,..,,,, NN, respectively. By comparing how indices are arranged, we can see that the en- tries of the matrix (No)) are essentially entries of the n th row of [Ej+alnjja inferring ( || S [[[,)+1),,,
As we have described before, Equation 6 is solved through multiple independent MMV problems. Thus, the error bound is also given regarding each individual MMV problem. Further, not only the jth step, but also the j + 1th step is solved through multiple independent MMV problems. Thus, we do not have the upper bound for Ej+1 in Equation 4 as whole but only the upper bound for each column block Ej+1Jn, 2" Fortunately, since all the noise entries in the jth step are contained as one single row of the noise block in the previous step, having the upper bound for each column block EEjt1]n,e is sufficient to derive the noise bound for the jth step, which is shown in Lemma 6 and Equation|
Now, we proceed to the proof of Theorem 2 Generally speaking, Theorem
2 is obtained by recur- sively applying Lemma 6 Particularly, focusing on the hierarchical sparse vectors, for the last step, i.e., the first mode unfolding, we solve
To = H1 + No))
leading I7 N2 SMV problems. They are SMV because there is only one column in each column block, and hence the MMV problem reduces the SMV problem. Lemma Commallinde iidiicates that
1 - UUlle < Σ ||ú,,, - [Ullnall & Σ
ofll(Ulbo2l + T EE l
All l1e + (a Uelns le + 72 Esla,lle)
Σ
ΣII (=1J=1
Fo U (0
II
E m besalle IlF
21

---

Published as a conference paper at ICLR 2026
We note that I + 1th level contains only block, leading [Er E₁. Using the relation U, U, + E, and Lemma leads || I <ail|Ur|le Tr ( |l.. This concludes the proof for s hierarchical sparsity. For standard sparsity, the upper bound for all standard sparse vectors the worst upper bound among all possible a corresponding to the sparsity level 8. Therefore, taking the maximum over Va E [s) concludes the proof.
For Kronecker-supported sparsity, since the support is shared among different blocks in the same level, it is unnecessary to introduce multiple MMV problems, but solve only one MMV problem. Icads to the final result. For the last step, we solve
Thus, recursively applying Lemma 5
To (1) + N(y:
which leads to the following relations,
1 - U. & a U Il + T N
s a ( + Ti (a 72 Neille)
a U.
7. N|
Thus, the proof is complete.
H ADDITIONAL NUMERICAL EVALUATIONS
H.1 COMPREHENSIVE STRUCTURED SPARSE VECTOR RECOVERY PERFORMANCE
This section presents a more comprehensive evaluation of our MSR framework compared to the state-of-the-art, consisting of complete results of Bection[]and new set of results where we vary the number of measurements with fixed SNR. We also include new metric named support recovery rate (SRR) defined as
SRR supp(()) supp(x:)
supp(#) supp(=)]]
where supp(.) returns the set of positions of the nonzero entries of the argument vector, | the cardinality of the argument set, & is the estimated sparse vector, and the ground truth.
returns
Table 4: Average runtime. A complete version of Table 1
Bold: the best result.
SNR
3 dB
dB
11 dB 15 dB
19 JB
23 JB
Recovery
sparse vectors
MSOMP-Seq
0.4256
0,4119
0,3827
0,3329
0.2204
0.0568
KroUMP
130.5405
108.0526
76.6942
39.9844
11.3774
0.7525
MSSBL-Seq
[.891
1.1016
0.5758
0.2218
0.1417
0.114]
MSSBL-P1
0.4517
4.2633
0,2658
1,8292
0.1531
0.1281
Recovery of a hierarchically sparse vectors
MSSBL-Seg
2.4900
2.0134
1,2501
0,1112
MSSBL-P1
0.4962
0.4607
1,3664
0.2513
0.1447
0.108]
MSHIP-9eq
0.0379
0,0305
0.0297
0.0247
0.0186
0.0168
HIGHTP
0.6512
0.5493
0.5204
0.5444
0.4398
0.4574
LTTP
2,2436
L,7770
1,3256
0.8450
0.8264
0.5311
MSHHT-Seg
D.0500
0.0510
0.0532
0.0509
0.0450
0.0434
IHI
8.2437
8.2412
8.2554
8.2917
8.2889
8.2789
Recovery
Kronecker supporrie sparse vectors
MSOMP
0.0042
0.0041
0.0040
0.0058
0.0026
0.00155
MSHTP
0.0011
0.0010
0.0010
MSSBL
0.0728
0.0587
0.0447
0.0279
0.0119
0.0051
SVD KrOSBL
37.1233
26.9816
14,2405
8.6036
5.4067
4.0681
AM-KrOBBL
55,9532
63.4676
75.9727
74.5840
51.7089
34,1331
HTP
0.9772
033377
0.4709
0.3466
0.2339
0.2323
MSIHT
0.0008
0.0007
0.0007
0.0007
0.0007
0.0007
IHT
6.0771
6.0760
6.0763
6.0690
6.0677
6,0535
KSITTP
0.1018
0.0730
0.0665
0.0881
We show a complete version of Figure 2 in Figure 7 We use Tensorlab (Vervliet et al. 2016) for tensor operation and Seq and P1. to represent the sequential and parallel (parfor function in
22

---

Published as a conference paper at ICLR 2026
105
10
0.5
1
1
10
10
SNR (dB)
SNR (dB)
(a) standard sparsity.
McSBL-Sa MESELP.
XSHTT--5
HIIP
RTP MSSHT-Saa IHT Sparsly -.ry/ Spansil
0.8
10"
NSE
10
0.6
10
10
10
20
25
10
20
25
SNR (dB)
SHR (dB)
(b) s hierarchically sparsity.
MEHIP MBSBL
ANG&nBBL
SVD-Kw6BL
v KSHTP Speaie
10a
10
0.6
100
10
10
(dB)
(dB)
(c) Kronecker-supported sparsity.
Figure NSE and SRR as functions of SNR. A complete version of Figure
Matlab [nc.. 2024)) implementation of Equation 5 they have the same recovery performance but different runtimes.
For fairness, we cap the number of EM iterations for SBL-based methods (MSSBL, AM-KroSBL, and SVD-KroSBL) two hundred, for HTP based methods one hundred, and for all IHT based methods to two hundred. For HTP- algorithms, we stop the iterations if the detected support remains the same two consecutive iterations (Foucari 2011). while IHT-based algorithms are ter- minated when the normalized difference between two consecuivv estimations is smaller than 10-6 For OMP-based algorithms, we stop when the norm the residual smaller than e/|y/|2- Here, the coefficient 0.05 is fixed for all OMP-based algorithms, which is empirically determined. We also prune small entries in hyperparameters for faster convergence for SBL-based algorithms.
In the recovery of s standard sparse vectors, compared to Figure 2a Figure 7a includes both the se- quential and parallel implementation of our MSSBL. Regardless of different runtimes as in Table 4 sequential and parallel implementations provide identical NSE and SRR results. Regarding runtime, MSSBL-P1 1 is only faster than MSSBL-Seq in low SNR cases. This is because in high SNR cases, the parallel overhead dominates, including data transfer and communication cost. As see in Fig-
23

---

Published as a conference paper at ICLR 2026
10
0.6
SAR
10
30
Ratio
Ratio
(a) 8 standard sparsity.
SSHET
x murr
A uSur
a1 9/ Spondty - Spanity
0.9
NSE
10
10
0.5
0.4
30'
0.35
0.4
0.45
0.5
0.55
0.65
0.7
0.75
0.8
0.35
0.4
0.45
0.5
0.55
0.6
0.65
0.7
0.75
(b) hierarchically sparsity.
A
Sponity Spuntry
10
0.6
SAR
30
10-4
0.35
0.55
0.65
0.7
0.35
0.4
0.45
0.5
0.55
0.65
Ratio
Rabio
(c)
Kromecker-suppooteed sparsity.
Figure 8: NSE and SRR as functions of the number of measurements.
24

---

Published as a conference paper at ICLR 2026
Table 5: Average runtime in seconds. Bold: the best result.
Recovery of spurse vectors
60
64
MSOMP-3ec
0,0735
0,0951
0.1268
0.1374
0,1937
0.2312
0.2918
KnOMP
0.7139
11941
21137
3,7952
7,4092
12.5439
22.4149
MSSBL-Seq
0.1726
0.1439
0.1364
0.1332
0.1368
0.1321
MSSBL-P1
0.1412
0.1464
0.1508
0.1517
0,1498
0.1527
0.3298
Recovery of hierarchically sparse vectors
M
48
52
56
60
64
68
72
MSSBL-Sog
0,2978
0.2594
0.2031
0.1692
0.1552
0.1407
0.1274
MSSBL-P1
0,1734
0.1614
0.1467
0.1363
0.1299
0,1248
0.1197
MSH
0,0204
0.0198
0.0190
0.0191
0.0178
0.0170
0.0168
HIRIP
0.3661
0.3622
0.3691
0.4189
0.4181
0,4599
0.5127
0.3071
0.3641
0.4733
0.4980
0.6407
0.6867
D
MSSIfF-seg
0,0527
0.0502
0.0479
0.0457
0.0449
0.0421
0.0416
TNT
4.7097
5,4549
6.2188
7.1292
8.3298
9.1231
10.1874
Recovery Kroneckor-supporid sparse vectors
M
12
13
14
15
16
MSOMP
0.0017
0.0017
0.0019
0.0020
0.0023
MSHIP
0.0010
0.0010
0,0011
0.0011
0,0013
MSSBL
D.0187
0.0141
0.0118
0.0092
0.0083
SVD-KroSBL
5,3134
4.8928
4.7247
5,0780
5,3094
AM-KroSBL
52,0891
5244799
49.6836
47.0243
42,6749
HTP
0.0715
0.1113
0.1658
0,2508
0.2745
MSIHT
0.0009
0.0008
0,0009
0.0007
0.0007
3,4581
4.0502
4,9489
6.0187
7.2859
KSHTP
0.1184
0.0908
0,1048
0,0885
0.0829
when the computation cost dominates, there is a significant gain in computation time, as a trade-off for memory usage. In the recovery of hierarchical sparse vectors, compared to Figure 2b Figure (alure/dese the performance of MSSBL. MSSBL exhibits a worse performance in low SNR scenario because it does not require the true sparsity level an input, while for IHT//TTP-based algorithms, this prior knowledge is necessary. However, MSSBL is still able to offer a comparable performance in high SNR scenarios, making it a powerful candidate when the prior knowledge is absent. In the recovery of s Kronecker-supported sparse vectors, compared to Figure we include IHT/HTP-basse algorithms in Figure 7c KSHTP is the algorithm we explained in Equation HTP Although the thresholding operator for Kronecker support is not optimal, KSHTP still offers the best SRR performance, followed by MSHTP. MSIHT has the least runtime, which is four orders less than its classic counterpart IHT. Overall, Figure and Table demonstrate that our MSR framework can offer similar or better performance with significantly reduced runtime.
ure
We next evaluate the performance of different algorithms by fixing the SNR and varying the number of measurements. The setting is as follows. For the s standard sparsity, we opt for H = _ with 1 = 2, and set M {48, 52, 5.. 72] and N = 80. The entries of H, and the nonzero entries æ are drawn independently from the standard normal distribution. We set s = 15, and the support is randomly drawn from uniform distribution. For hierarchically sparse vectors, also opt for = 2, and set M = (48, 52, 56, ... .72). = 80, and s = 15. the Kronecker- supported sparsity model, we opt for 1 3, and set M {12, 13,.... 16). N 18, and 4. We adopt the additive white Gaussian noise with zero mean with SNR (dB) = 20. Ratio is defined as M/Ñ -I' Mj/N, (M/N)'. We consider NSE, SRR, and runtime for performance evaluation. We follow the same way to cap the number of iterations. Results in Figure Table are obtained through two hundred independent trials. Overall, we observe similar trends as in Figure and Table Our MSR is able to provide comparable or better performance with reduced runtime, demonstrating the efficacy of exploiting the Kronecker product structure in the recovery process.
H.2 COMPARISONS WITH TRADITIONAL COMPRESSED SENSING APPROACHES
In this section, compare our MSR to traditional compressed sensing algorithms, including IHT. HTP, SBL, OMP, and the & norm-based basis pursuit denoising (BPDN) (Foucart & Rauhut] 2013) (basisPursufs function in Matlab (Inc.] )024
25

---

Published as a conference paper at ICLR 2026
4 BSUP
Amm
KHOW pair QUTE
q COUP
Oreentl 4-300 QJECITY DAY Q ARTT
ANn
XNCMMP Hoir to Quss gini +usin pir! qUNI'
16
10
19
10"
10
10
(di))
(si))
(ci))
(a) Standard sparsity.
(b) Hierarchical sparsity.
(c) Krooecker-supportee sparsity.
Figure 9: NSE as function of SNR compared to traditional compressed sensing algorithms.
For all three models, we set M, M, N, = N, and = for i E [7]. For the standard sparsity, we opt for H , H2 with T= 3, M 12, and N 15. The entries of H, and the nonzero entries of
drawn independently from the standard normal distribution. We set = 8, and the support is randomly drawn from a uniform distribution. For & hicrarchically sparse vectors, we opt for I 2. MM 35, N 40, and s 8 per dimension. Here, supports generated by first selecting blocks uniformly at random, then assigning support within each block uniformly. In the Kronecker-supported sparsity model, models, the measurement noise is zero mean white Gaussian noise whose variance is determined of {0, 5, 10, 15, 20, 25}. We follow the same condition to cap the iterative algorithms, and for BPDN, we cap the number of iterations at fifty. Compared to Section | and Appendix H.T we downsize the measurement matrices mainly for computational feasibility. With the same condition in Appendix H.T it is hard to evaluate traditional algorithms such as BPDN and SBL. The NSE shown Figure 9 are median and 25%/75% quartiles, while Table shows the average runtime, both over fifty independent trials.
opt for I
M 12, N 15, and 4. In all
by SNR (dB)= 10 1og10 [E{ H a /// // [[]
As shown in Figure 9a HTP achieves the lowest runtime and high accuracy but relies on prior knowledge of the true sparsity level 8, which is generally unavailable in practice. While SBL yields higher reconstruction accuracy than MSSBL in the high SNR regime, it incurs a runtime two or- ders of magnitude higher, limiting its scalability. MSSBL emerges the most robust solution: it outperforms traditional methods (IHT, OMP, BPDN) and has lower runtime than SBL., offering a balance between reconstruction accuracy and computational efficiency, without requiring specific prior knowledge.
We illustrate the hierarchical sparsity recovery against traditional compressed sensing algorithms in Figure 9b] Compared to HTP, MSHTP achieves almost the same reconstruction accuracy but with to two orders of less runtime. This also happens MSSBL compared to SBI. MSSBL, in this case, is also a balanced option when the true sparsity level is unknown, without sacrificing efficiency significantly. OMP and IHT are only slightly worse than their counterparts, i.c., MSOMP and MSIHT, but the gain in runtime is significant by two orders of magnitude. Finally, Figure 9c contains the results for -upppoor- sparsity recovery against traditional benchmarks. MSSBL constantly achieves the best performance, seconded by MSHTP and MSOMP (high SNR case), with two to three orders less runtime than their counterparts.
H.3 SPARSE VECTOR RECOVERY PERFORMANCE WITH VARYING NUMBER OF DIMENSIONS I
MASBBL-Sa x X. DASSHT--S
MSSTT- PWA.SIIIITT--S-S- AASIITT-P
10
10
101
10
NBF
10
10
10
10
2
3
2
3
System Order
System Ordes
Figure 10: NSE and runtime of MSR functions of system order I.
26

---

Published as a conference paper at ICLR 2026
Table Average runtime for comparison with traditional compressed sensing algorithms. Bold: the best result.
d
10
dB
Recovery of sparse voctors
MSOMP
0.0607
0.0578
0.0095
MSSBL
0.4547
0.4096
0.3356
0.1929
0.0969
0.0583
HTP
0,0210
0.0129
0.0100
0.0098
0.0091
0.0109
III
2,2580
2,2627
2,2629
2.2577
2,2473
2.2425
OMP
22.4379
20.1685
15.7469
8,9826
2,3411
0.0363
SBL
37.9208
21.5275
13,2480
10.9864
10.0749
9,4921
BPDN
102.8134
102.0363
T01.4731
99.6498
96.1337
85.2775
Recovery hierarchically sparse
MSOMP
0.0246
0.0248
0.0232
0.0233
0.0123
0.0024
MSHIP
0.00447
0.0037
0.0031
0.0032
0.0033
0.0028
MSIHT
0.0066
0,077
0.0082
0.0088
0,0889
MSSBL
0.4694
0.3997
0.2533
0.1538
0.0424
0.0172
HTP
0. 0969
0.1225
0.0777
0.0620
0.0363
0,0274
0,4834
OMP
6.5255
5,7934
4,4745
2,6242
0.7864
0,0365
SBL
18,3118
18.1878
8.3721
3,2350
2.1724
1,8597
BPDN
20.9009
20.8322
20,3150
19.2756
17,9916
16.0579
Recovery of Kronecker-uupooee sparse vectors
MSOMP
6.0031
0.0027
0.0027
0.0025
0.0023
0.0014
MSUIT
0.0011
0,0010
MSIHT
0.0010
0.0009
0,0010
0.0008
0,0010
0,0009
MSSBL
0.0540
0.0500
0.0511
0.0322
0.0183
0.0057
HIP
0.3165
0.2945
0,2105
0.0945
0.0678
0.0344
THT
1,8995
1,9024
1,8995
1,899
1.152
1,8953
OMP
23.6829
21.5929
16.4708
9.4103
2.68339
0.1793
SBL.
46.1748
42.6398
29.0591
15.0069
11.0725
9.5650
BPDN
106.9569
Here present results where fix the size of each factor matrix H, but vary I from 1 to 4, to demonstrate the ability of our MSR to handle arbitrary system order I.
We consider N₂ N and M, M for i [4], opt for N 50, and determine M and through [0.4N]] as in Section We fix SNR at 20dB and results in Figure 10 shows how the NSE and runtime (median with 25% 25% quartiles) of MSSBL, MSHTP, and MSIHT scale with different for Kronecker-supported sparsity model.
M [(0.6N) and
In Figure TO we observe similar trends as in Figure 3 Parallel and sequential implementations have identical recovery performance but different runtimes. Thus, for NSE, we only show the results for parallel implementation. The only exception in Figure I0 compared to Figure 3 is that for = 2, parallel implementation for MSHTP and MSIHT requires more computation time than the sequential implementation. This is because when the problem is computationally light, parallel overhead dominates the time consumption, including data transfer and communication, rather than computation itself.
H.4 APPLICATION: CHANNEL ESTIMATION FOR INTELLIGENT REFLECTING CRRFACE-AIDER WIRELESS SYSTEM
Channel Estimation NSE
10
SNR (dB)
AM-KreSBL
SVD-KroSBL
MSSBL
5
138.9725
3.1272
0.0212
10
121.7270
2.8479
0.0147
10-4
15
77.7143
2,7911
0.0099
AM KreSBL MISTEEL
20
43.9860
2,7752
0.0066
25
26.3338
2,7479
0.0052
10
SNR (dB)
18.7202
2.7450
0.0052
Figure 11: NSE of different schemes.
Table 7: Average runtime in seconds.
27

---

Published as a conference paper at ICLR 2026
An intelligent reflecting surface (IRS) is a reconfigurable meta-surface consisting of large num- ber of adjustable reflecting elements. By changing the reflection coefficients of elements, IRS can reflect the signal to certain area to improve the coverage of the wireless communication systems frequency bands. But properly configuring the IRS requires the channel state information. It is hence important to develop efficient and accurate channel esti- mation algorithms for the IRS-aided system. In this section, we consider an uplink narrowband IRS- aided MIMO system, consisting of a I walf-wavelength spacing-antenna transmitter mobile station (MS), an R balf-wwaelength spacing-element uniform linear array IRS. If the channel matrices of the MS-IRS, and the IRS-BS channel are denoted as 4sss E C and <p E CRx1, respectively, according to the geometric channel model (You et al.] 2022
operating at pallimeteerwave/trrahhe
-antennnannnnnaaaaaa receiver base station (BS), and an L half- wavelength
Wang et al. 2020) Alkhateeb al., 2014 He & Joseph 2023
PsS E CLXT and 4Bs E
CAXL
can be formu lated as
PME
LT Pus
BMS y (MMs, (axs))
(8)
Pas
RL
PB
BBS (a(Bs,
(9)
PBS
where PMs and Pas are the number of paths between MS and IRS, and IRS and BS, respectively. The angles ΦMS.p- OMS, OBS..p, and øBs represent the pth AoA of the IRS, and the pth AoD of the MS. the pth AoA of the BS, and the AoD of the IRS, respectively, while Bys.. and $Bs,p are the complex path gains. Steering vector (1) E CO for any integer Q and angle & is defined e cos T1T. Then, the cascaded MS-IRS-
as (v) L/VQ[L, co given by ps diag(( ( Ms for a given IRS configuration 0 C² whose ith entry of models the reflection of the ith IRS element. Channel estimation problem targets to estimate the cascaded chan- nel Bs diag ( 4 MS given any 0, which is sufficient for subsequent tasks such as beamforming
cos
channel is
Waan et al. 2020).
Wireless channels operating on millimeter-waae//rrrheeet bands are intrinsically sparse due se- vere path loss. To reveal this sparsity, we adopt three sparsifying bases An. AL, and Ar, corre- sponding to the angular domain of the array at BS, IRS, and MS, respectively. Such bases contain steering vectors evaluated over N grid angles (1², IN n=1 such that cos((),) = 2n/N 1 MMa et al., 2022), defined as Ag = [aq(v)), (12))
ag [WN]] E CQXN for any integer Q > 0. Then,
Equation 8 and I reduce to
BS AR AH and 4MS- = AI AH
(10)
where vectors FR. LL,dd LL,a; TT E CN are the unknown channel representations over the known sparsifying bases. They are sparse due to the intrinsic sparsity of the channel.
Channel estimation is performed by processing the received pilot signals, given the knowledge of the sent pilot signals and the training IRS configurations. Suppose we allocate time slots for channel estimation, over which the channel is considered to be constant. We vary IRS configurations for Ki times, and for each different configuration, we transmit the same set of pilot signal G E CTxKp over Kp time slots such that K = KjKp. The received signal Ys II Kp corresponding to the kthh training configuration 02 is
Y Dns diag(Ox)) P GNN
(11)
where N& E CI is the noise. Substituting Equation 10 into Equation 11 and vectorizing the received signal { followed by some algebraic operations He & Joseph 2023), we have
y (Hr @ HT 8 HR) x + n Hæ + n E CKK
(12)
where H. E C is formed by the first N columns of ( Φ A whose N2 columns are just N repetitions the columns of H., H = AG. and HR = AR. with @ being the Khatri- Rao product. We collect K1 IRS configurations {0x] in matrix O E Cl and define 31. @ wi @ ØR E C with कt. E C being the scaled version of the first N entries of Ft.a L corresponding to the removal of redundant columns in θ T(AA] @ AH)T. Denoting
28

---

Published as a conference paper at ICLR 2026
H1 H3 with Kj as My. HT H2 with Kp M2. HR H1 with R as M1, and N = N₂ N2 = Ni, the channel estimation problem is transformed into a KCS problem following the form of Equation I with unknown
($3, $2 $1) (PMs: 1, PBS) Kronecker-siv suppereerr-upprteed sparse vector æ. With the estimated &; the estimated channel for given IRS configuration O, i.e. Bss diag(0) @ MS. obtained by reshaping (& Ms @ $Bs) / with known size R x T, where + I O Bs with known size RT X L is reconstructed as vec( rec( MS 3 PBB) A @ A Ø AR) ) with & being the first N columns of [A]. @ A
We set T == R 16, and I, 256. To set up the sparsifying bases, opt for N 18. Pilot signals contained in G are randomly generated quadrature phase shift keying symbols using a uni- form distribution while the IRS training configuration randomly drawn from uniform distribution (+1/VL) <Lin et al. 2021). Regarding the pilot signals and training IRS configurations, we con- sider K1 = -p I0, making H3 E C10x18, 18, H₂ E C10x 18, and H1 E C C16x 18 To model the scatters, we set PMS PBS 3 and all angles OMS.p. aMSs QBSS,p. and ÓBS are drawn uniformly and independently from the grid points, while path gains fMs. and $BS,, are drawn independently from complex standard normal distribution (He & Joseph 2023). We compare MSSBL. with SVD- [AM-KroSBL in He & Joseph (2023) with the same way to cap iterative algorithms. Metries include channel estimation NSE given by K, EXE [ ding(0, ) ass diag(O, (PMss and runtime.
|| p diag((0, s
Figure 11 shows 25%/50%/75% quartiles of NSE, while Table 7 shows the average runtime, both over fifty independent trials. We observe that all three algorithms provide comparable channel esti- mation performance, while MSSBL has two orders less runtime than EVD-KroSBL. and four orders less runtime than AM-KroSBL, making it more efficient in this application scenario.
H.S APPLICATION: FOREMAN VIDEO SEQUENCE RECOVERY
35
Traditional
MSR
X Traditioall CS OSn
0.05
77.5294
10.5590
30
0.10
91.2761
9.0006
0.15
70.1564
9.7536
PSNR
0.20
75.8293
9.9201
0,25
67.9980
11.6994
0.30
99.2622
11,1384
20
79.0998
IIIIII
0.40
109 8495
11.1645
15
0.3
0.45
92.8983
7.7705
Ratio
66.3104
6.7719
Figure 12: PSNR of different schemes.
Table 8: Runtime in seconds.
This section presents the recovery results for Foreman video sequence to demonstrate MSR's ability to deal with real-world dataset and its superiority over the traditional compressed sensing algo- rithm. In the experiment, we follow the settings in Duarte & Baraniuk (2011a and to make this paper self-contained, we provide a brief overview of these settings.
Tested frames are generated by cropping around the center to form frame size of 128 x 128 pixels and there are in total eight frames used in the experiment. To spatially sparsify the image content within a single frame, we vectorize each frame and adopt a 2D inverse discrete wavelet transform basis W1 E R 16384x16384 applied to the sparse coefficients of each video frame. For the sparsity in the temporal dimension, we turn to an 1D inverse discrete wavelet transform basis W2 E 8 It exploits the correlation between frames to sparsify the signal over time. Suppose the video sequence is denoted by 0 E R 131072 where 131072 = 8 (frames) x 128 (row pixels) x 128 (columns pixels), then its relation to sparsifying bases and the sparse coefficient vector x E R 131072 is 0 (W2 @ W1).. To compress the vidco sequence, we x that there is no temporal compression while S E M 16384 from a subsampled permuted Hadamard transform denotes the spatial compression. Here, M is the number of measurements taken in one frame. This leads to the following meassurement/sparsffyinng model
mcasurement matrix as I 2 S, where I means
y (W2 @ SW1)..
(13)
The goal is to obtain sparse coefficient x using the compressed measurement y. and finally we reconstruct the video sequence as 0 (W₂ Wi
29

---

Published as a conference paper at ICLR 2026
Figure 13: Comparison of reconstructed Foreman video frames. Row 1: ground truth. Row 2: traditional CS with 0.1. Row 3: MSR with Ç 0.1, Row traditional CS with Ç 0.5. Row 5: MSR with Ç - 0.5.
Denoting W2 as H2, SW1 as H1, and 2, we note that Equation 13 is mathematically equivalent to Equation hence be solved using our MSR. The benchmark in this experiment is & based basis pursuit (Duarte & Baraniuk 201 where Equation [3 is treated as traditional compressed sensing problem and the Kronecker structure of W2 @ SW, is ignored. For fair comparison, we adopt the same & solver (van den Berg & Friedlander 2008 2019) in MSR as Duarte & Baraniuk (2)11a) with the same stopping criterion. We determine the number measurements M as M- [6384ç] with measurement ratio Ç E {0.05, 0.10, ... .,445, 0.50}. We use the peak signal-to-noise ratio (PSNR) and runtime evaluation metries. Results are shown in Figure [2] and Table 8] where traditional CS refers to traditional compressed sensing algorithm. We observe that MSR achieves the same PSNR as that traditional CS with roughly one order magnitude less runtime, effectively demonstrating the efficacy of our MSR on real-world dataset. We also compare all used video frames in Figure 3
30


Reviews:
[{"review_id": "AXIEjnRrxK", "reviewer": "Reviewer_27dS", "summary": "The paper studies the performance of Kronecker compressed sensing (CS) for three different hierarchical/structured sparsity patterns. New algorithms that mix reshaping of signal, observation, and sensing tensors and joint sparsity algorithms applied to tensor slices allow for improvements in performance versus agnostic baselines.", "questions": "Are there applications where the newly proposed hierarchical sparsity structures arise naturally?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "The paper covers theoretical analysis, algorithmic contributions, and numerical performance for several different signal setups. \n\nThe analytical results expand over existing contributions in Kronecker CS by proposing and leveraging custom restricted isometry properties, obtaining better scaling laws for the number of measurements necessary from random CS matrix constructions. The results also include the commonly observed tolerance to the presence of noise in the measurements.\n\nThe proposed algorithms exploit the structure of Kronecker matrices for CS and sparsity to reduce the computational complexity of recovery vs. their agnostic counterparts.\n\nThe numerical comparisons are convincing.", "weaknesses": "It is not always clear if these proposed sparsity structures find readily available applications where Kronecker CS is feasible and an improvement over baselines. All experiments are based on synthetic data.\n\nThere appears to be no discussion of the role of the sparsity transformation, which implicitly would have to follow a Kronecker product formulation as well.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "p8Qzx5SzTx", "reviewer": "Reviewer_64D7", "summary": "This paper proposes a multi-stage sparse recovery framework for Kronecker compressed sensing, tailored to three sparsity models, with theoretical guarantees and simulations showing comparable performance and faster runtime than state-of-the-art methods.", "questions": "Besides the two major concerns mentioned above, I have the following questions:\n\n1) The paper introduces three sparsity models, but it does not specify which types of real-world data these models are tailored for. Despite asserting broad applicability, the paper lacks simulations using authentic datasets. \n\n2) As previously mentioned, a critical comparison with traditional compressed sensing methods in terms of reconstruction accuracy is essential, given the high computational complexity of KCS.\n\n3) In the experiments, the analysis is limited to cases with a number of stages no greater than three, namely $I\\leq 3$. Why was this restriction imposed? Ideally, the performance variations across different stage numbers should have been explored.", "limitations": "", "rating": 4, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "1）Developed a versatile multi-stage sparse recovery algorithmic framework tailored to three different sparsity models (standard, hierarchical, and Kronecker-supported) for the Kronecker compressed sensing problem.\n\n2）Provided theoretical recovery guarantees by analyzing the restricted isometry property of Kronecker product matrices under different sparsity models, and demonstrated comparable recovery performance with significantly reduced run time through simulations.", "weaknesses": "1）Compared to traditional compressed sensing, Kronecker Compressed Sensing (KCS) has not garnered sufficient attention, primarily due to its excessive complexity—including challenges in measurement matrix construction, signal reconstruction and hardware implementation—while lacking significant performance advantages that could justify its complexity. One notable strength of the method proposed in this paper is its reduced reconstruction complexity. However, we are particularly interested in how it compares to traditional methods in terms of runtime and reconstruction accuracy, yet the paper appears to lack relevant experimental evaluations.\n\n2）The construction of measurement matrices remains a major challenge in KCS. While the paper claims to improve Restricted Isometry Property (RIP) conditions, there seems no  concrete method for constructing matrices that satisfy these enhanced conditions.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Wg8kHeQSQO", "reviewer": "Reviewer_Wv5h", "summary": "The authors consider the problem of recovering sparse vectors $x$ from Kronecker measurements, i.e., measurement matrix of the form $H=H_1\\otimes \\dots \\otimes H_N$. The overarching idea is to perform the recovery in stages: By unfolding of the measurement $y=Hx$ in different manners, the recovery problem can be viewed as a number of  MMV problems of the form $y= H_j U_j$, that can be solved in an iterative fashion. They coin this class of algorithms MSR, multi-stage recovery. Importantly, their method is improved if the $x$ also is assumed to have a Kronecker sparse structure, since the MMV problems become coupled.\n\nThey derive a convergence guarantee for their method that is valid under the assumption that each $H_i$ has the RIP. Their bound can be applied to different types of sparsity in $x$, and in particular improves when $x$ Kronecker-sparse. It however potentially suffers from the potential of error propagation between the folded stages.", "questions": "** Relation of RIP analysis to [Roth et al; 2020]** Although the formulation of the results are different, it seems like the formulation and proof of Theorem 1 in the manuscript at hand is essentially the same as the formulation and proof of Theorem 4. Can the authors comment on this, and in particular pinpoint exactly how the analyses are different?\n\n** Complexity discussion and HiHTP ** In the complexity discussion, the authors claim that the space- and time-complexity for HiHTP ($I=2$) both are $M^2N^2$. However [Roth et al; 2020] claim in their paper that the per-iteration time- complexity of their algorithm (for generic measurement matrices) is (in this notation) $MN^2$. Can the authors clarify?", "limitations": "", "rating": 6, "confidence": 4, "soundness": 2, "presentation": 2, "contribution": 3, "strengths": "The main strength of the paper is that it proposes a method which performs better for Kronecker-based sparsity than for general hierarchical sparsity, which with the exception of He and Joseph, 2025, seems to be novel. \n\nThe recovery guarantee that the authors provide is interesting in that it only requires each constituent matrix to have a constant RIP ($\\delta_{s_i}(H_i)<1/\\sqrt{3}$) for convergence in the noiseless case. This is actually a lot better than previous guarantees -- using the bound in Theorem 1 to get an RIC (which is what is needed for e.g. the HiHTP to converge), one essentially need the sums of the $\\delta_{s_i}(H_i)$ to be less than $1/\\sqrt{3}$, which intuitively means that each constant must be smaller than $1/\\sqrt{I\\cdot 3}$, which will be very low for a high number of Kronecker levels. It is also interesting to see that their technique, through Corollary 1, improves on the previously quite weak bounds for standard sparsity.", "weaknesses": "The paper is sometimes hard to follow, mainly due to the complexity of the notation related to tensor-algebra and the different sparsity notions. I think that refinement of the illustrations would go a long way to alleviate the problems. The color coding in Figure 4 is in my opiniion rather confusing than helpful.\n\nThe RIC analysis offers some new interesting perspectives, but it seems very close to the one performed in [Roth et. al; 2020] and [Roth et al; 2018]. The relations between the two should be better explained, see questions below.\n\nThe superiority of the proposed methods compared to the SOTA, both for hierarchically sparse vectors and Kronecker-sparse vectors, are mainly argued for via comparison of run-times. This metric depends a lot of the specifics of the implementation. Could something more objective be measured and presented, such as the number of iterations? I particularly stumble upon this since the discussion regarding time-complexity on page 5 seems to contain mistakes for HiHTP, see below.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "5d5fj6GNgH", "reviewer": "Reviewer_MA64", "summary": "In the Kronecker compressed sensing (KCS) framework, the measurement matrix is composed by the Kronecker product of multiple factor matrices as $ H = H_I \\otimes \\cdots \\otimes H_1$. While this allows for efficient measurement of signals with multidimensional structures, the dimensionality increases rapidly (e.g., $O(N^I))$, which makes existing methods computationally expensive. They also fail to fully utilize the individual factor matrix structures and hierarchical sparsity.\nThis paper addresses these difficulties by introducing the {\\em hierarchical view}. \nThe contribution of this paper is composed of the following three parts. \nFirst, it is theoretically demonstrated that the Kronecker product structure is a \"hierarchical block partition\" and that each factor matrix $H_i$ measures the sparsity of a different \"hierarchical level.\" Taking advantage of the Kronecker structure, an algorithm, Multi-Stage Recovery (MSR), that performs reconstruction sequentially (or in parallel) for each layer is proposed. This reduces the computational cost for signal reconstruction to $O((MN)^I)$ to $O(MN^I)$. \nSecond, a theoretical guarantee is developed by introducing a new type RIP-based condition. \nFinally, for all three sparse models, MSR methods (MSOMP, MSHTP, MSSBL, etc.) show comparable or superior accuracy to existing methods (KroOMP, HiHTP, KroSBL, etc.), while reducing computational time by 1–3 orders of magnitude.", "questions": "This paper refers to radar imaging and wireless communication as examples of applications where the proposed method is useful, because the sparsity of each signal can be separated into dimensions. Can the usefulness of this method be verified for real data in these fields?", "limitations": "", "rating": 8, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "By introducing a new perspective into the Kronecker compressed sensing framework, we propose a signal restoration algorithm that can significantly reduce the required computational effort compared to conventional methods. At the same time, we also conduct novel performance evaluations, theoretically guaranteeing performance, and demonstrate the practical usefulness of the proposed algorithm through numerical experiments.", "weaknesses": "Kronecker compressed sensing is only possible when the observed signal has multidimensional structural sparsity or separable statistical structure, which limits the practical applications in which the proposed method is useful.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)