Minimax estimation in linear models with
unknown finite alphabet design

arXiv:1711.04145v1 [] 11 Nov 2017

Merle Behr∗ and Axel Munk∗,†
Institute for Mathematical Stochastics∗ ,
University of Göttingen,
and
Max Planck Institute for Biophysical Chemistry† ,
Göttingen, Germany
Email: {behr,munk}@math.uni-goettingen.de
We provide minimax theory for joint estimation of F and ω in linear
models Y = F ω + Z where the parameter matrix ω and the design matrix
F are unknown but the latter takes values in a known finite set. This
allows to separate F and ω, a task which is not doable, in general. We
obtain in the noiseless case, i.e., Z = 0, stable recovery of F and ω from
the linear model. Based on this, we show for Gaussian error matrix Z that
the LSE attains minimax rates for the prediction error for F ω. Notably,
these are exponential in the dimension of one component of Y . The finite
alphabet allows estimation of F and ω itself and it is shown that the LSE
achieves the minimax rate. As computation of the LSE is not feasible, an
efficient algorithm is proposed. Simulations suggest that this approximates
the LSE well.
Keywords: Finite alphabet, Combinatorial linear model, Lloyd’s algorithm, Least
squares estimator, Minimax estimation, Exact recovery, Blind source separation, Measurement uncertainty.
MSC 2010 subject classifications: Primary 62F12, 62H30, Secondary 62F30, 62J05.

1 Introduction
Linear models (LM) are arguably one of the most prominent tools in statistical modeling. In such a (multivariate) LM one observes a matrix Y ∈ Rn×M ,
Y = F ω + Z,

(1)

which is linked to the parameter matrix of interest ω ∈ Rm×M via a design matrix
F ∈ Rn×m and an additive noise matrix Z ∈ Rn×M , which in this work is assumed
to be i.i.d. Gaussian Zij ∼ N (0, σ 2 ) for i = 1, . . . n, j = 1, . . . , M , with (unknown)
variance σ 2 > 0. Usually, the design matrix F is assumed to be known and analysis
on ω is performed conditioned on F , as, e.g., in classical ANOVA. In contrast, in the
following, we want to consider a specific LM where the design matrix is unknown and
has to be estimated from the data Y jointly with the parameter matrix ω. In general,
separation of F and ω from F ω is not possible, of course, and therefore existing

1

approaches (see, e.g., Pananjady et al. (2017)) focus on estimation of F ω. However, if
we assume that F can only attain values in a known, finite set A = {a1 , . . . , ak } ⊂ R,
denoted as finite alphabet, we will show that the LM in (1) becomes identifiable, that
is, F and ω can be separated from F ω, under rather weak assumptions on ω and F .
The aim of this paper is to provide estimates and algorithms for joint recovery of F
and ω and to develop statistical theory for these quantities. For better understanding,
it is convenient to rewrite the LM (1) with unknown design matrix F as a blind source
separation problem (the terminology is borrowed from the signal processing literature)
Y·l =

m
X

F·i ωil + Z·l ,

l = 1, . . . , M,

(2)

i=1

with m source signals F·1 , . . . , F·m ∈ An , each only taking values in the finite alphabet
A, which are obtained in M mixtures with unknown mixing weights ω·1 , . . . , ω·M ∈
Rm . Therefore, we will denote model (1) with F ∈ An×m unknown as the Multivariate finite Alphabet Blind Separation (MABS) model. Throughout the following bold
letters, e.g., F , ω, denote the underlying truth generating the observations Y in (1).
Further, throughout the following kxk denotes Frobenius norm (for matrices x) and
Euclidean norm (for vectors x).
Applications MABS occurs in many different fields. For instance, in digital communications (Proakis, 2007; Talwar et al., 1996; Verdu, 1998; Zhang and Kassam, 2001;
Sampath et al., 2001), where m digital signals (e.g., binary signal with A = {0, 1})
are modulated (e.g., with pulse amplitude modulation), transmitted through wireless
channels (each having different channel response), and received by M antennas. In
signal processing this is known as MIMO (multiple input multiple output) and (ignoring time shifts, i.e., considering instantaneous mixtures) can be described by MABS
when the channel response is unknown, see e.g. (Talwar et al., 1996; Love et al., 2008)
for details.
Another example where MABS is relevant arises in cancer genetics (Yau et al.,
2011; Carter et al., 2012; Liu et al., 2013; Ha et al., 2014). Specific mutations in
cancer tumors are copy number variations (CNV’s) where some parts of the genome
are either duplicated or deleted. CN’s of a single tumor only take integer values, i.e.,
A = {0, 1, 2, . . . , k} (with good biological knowledge of a maximum copy number k,
see e.g., (Liu et al., 2013)). However, tumors are known to be heterogeneous, i.e., they
consists of a few different types of tumor cells, so called clones, see e.g., (Shah et al.,
2012; Greaves and Maley, 2012). In whole genome sequencing (WGS) data the CN’s
of the single clones overlap according to the relative (unknown) proportion of the clone
in the tumor. Important for the analysis of this data is that often M different probes
of the tumor cells are available, taken at different time points or at different locations.
Each of these contain the same clones but at different relative proportions. This can
be modeled with MABS, where m is the number of clones, M is the number of probes,
ωij corresponds to the unknown proportion of clone i in probe j of the tumor, and
F·i corresponds to the CN’s of clone i, see e.g., (Beroukhim et al., 2010; Greaves and
Maley, 2012; Shah et al., 2012) for details. For M = 1, see (Behr et al., 2017). Analog,

2

one can model point mutations in tumors with MABS with alphabet A = {0, 1}, where
Fji = 1 if and only if point mutation j is present in clone i (Jiang et al., 2016).
Two simplifications Motivated from the application in cancer genetics, where the
mixing weights correspond to physical mixing proportions of DNA strands, we will
in the following assume that the mixing weights ωij are positive and sum up to one
for each j. This assumption simplifies the corresponding identifiability conditions to
decompose F and ω uniquely. However, we stress that all results can be extended for
general mixing weights which only requires a slight modification of the corresponding
identifiability assumptions, see Remark 2.4. The minimax rates which are derived in
the following do not depend on this assumption. More precisely, for a given number
of sources m and a given number of mixtures M , the set of possible mixing weights ω
is defined as
Ωm,M :=
(
ω∈

Rm×M
+

: 0 < kω1· k < . . . < kωm· k,

m
X

)
ωij = 1 ∀j = 1, ..., M

.

(3)

i=1

Note that a fixed ordering of the row-sums is necessary as otherwise for any permutation matrix P one finds that F ω = F P P −1 ω with ω and P −1 ω both valid mixing
weights. Moreover, throughout the following, we may assume w.l.o.g. that the fixed
given alphabet A is ordered and that a1 = 0 and a2 = 1, that is
A = {0, 1, a3 , . . . , ak }

with

1 < a3 < . . . < ak .

(4)

Otherwise, one may instead consider the observations (Yj −a1 )/(a2 −a1 ) with alphabet
−a1
−a1
, . . . , aak2 −a
} in (1).
A = {0, 1, aa23 −a
1
1
Identifiability A minimal requirement underlying any recovery algorithm of F and ω
from (a possibly noisy version of) G := F ω in (1) to be valid is identifiability, that is,
a unique decomposition of the mixture G into finite alphabet sources F and weights
ω.
For illustration, consider a binary alphabet A = {0, 1} with two sources m = 2 and a
single mixture M = 1. The question is as follows: Is it possible to uniquely recover the
underlying weights ω ∈ R2+ and sources F ∈ An×2 from the mixture G = F ω ∈ Rn ?
In this example the answer is simple: As the entries in F can only attain the values
0 and 1 the smallest possible value for Gj , j = 1, . . . , n, is 0, which corresponds to
both sources taking the smallest alphabet value Fj1 = Fj2 = 0. Analog, when Fj1 = 0
and Fj2 = 1, Gj takes the second smallest possible value, denoted as ω1 (recall that
0 < ω1 ≤ ω2 and ω1 + ω2 = 1 by (3)). Similar, the third smallest value for Gj then
equals ω2 with Fj1 = 0, Fj2 = 1 and the largest value equals 1 with Fj1 = Fj2 = 1.
Thus, one can (almost) always uniquely identify ω and F from G = F ω. There are
just two situations where this fails:
1. If ω1 = ω2 , one cannot distinguish from G whether Fj1 = 0, Fj2 = 1 or Fj1 =
1, Fj2 = 0.

3

2. If Fj1 = Fj2 , one cannot distinguish ω1 , ω2 .
Consequently, in order to guarantee identifiability (in this simple example), we have
to exclude these two situations. That is, we need to exclude from the parameter space
the single weight vector ω = (0.5, 0.5) (the only one with ω1 = ω2 ) and sources F =
(F·1 , F·2 ) with equal components F·1 = F·2 (or equivalently ω1 , ω2 6∈ {G1 , . . . , Gn }).
Clearly, this is not very restrictive in most situations, as it simply excludes that only
one source is visible.
Now we turn to the general case, of arbitrary A, m, and M . It is shown in (Behr and
Munk, 2017) that identifiability has a complete combinatorial characterization via the
given alphabet and that the above assumptions can be extended to a universal (for
any A, m, M ) simple sufficient condition, called separability, which then guarantees
identifiability. For sake of completeness and to ease reading, we recapitulate these
conditions in the following in our context.
First, we discuss conditions on ω. For fixed ω the mixture G = F ω can take
Pmany of
at most k m (recall that the alphabet A has size k) values of the form eω = i=1 ei ωi
with e = (e1 , . . . , em ) ∈ Am (elements in Am are considered as row vectors). Clearly,
if for any two e 6= e0 ∈ Am it holds that eω = e0 ω, then F is not identifiable from G,
in general, as it cannot be distinguished from G whether F = e or F = e0 . Hence, we
require that the alphabet separation boundary (ASB) (Behr and Munk, 2017), that is,
the minimal distrance between any of these values, is positive, i.e.,
ASB(ω) = ASB(ω, A) :=

min

e6=e0 ∈Am

keω − e0 ωk > 0.

(5)

Further, we will see that when eω is corrupted by noise as in the MABS model (1),
the ASB quantifies perturbation stability of G when any two of these values are close,
that is, keω − e0 ωk < δ for small δ > 0.
Recall that for Ωm,M in (3), to ensure identifiability one has to assume kω1· k < . . . <
kωm· k, as otherwise the ordering of the rows of ω (and columns of F , respectively)
are not well defined via the mixture F ω. Hence, in a noisy setting, as in the MABS
model, the minimal distance between the row norms of ω will be crucial for separation
and hence recovery of F and ω. Therefore, for ω ∈ Ωm,M we introduce the weights
separation boundary
W SB(ω) = W SB(ω, A) :=

1 + mak
2

min (kωi· k − kωi−1· k) .

i=2,...,m

Figure 1.1 illustrates for m = M = 2 in (3) the set



a
b
Ω2,2 =
: a, b ∈ [0, 1), a + b < 1 .
1−a 1−b

(6)

(7)

For A = {0, 1}, regions where min(ASB(ω), W SB(ω)) ≥ δ are displayed blue in
Figure 1.1, for different values of δ. In particular, this illustrates that the set of
regularized weights {ω ∈ Ωm,M : min (ASB(ω), W SB(ω)) ≥ δ} is reasonably large
when δ is sufficiently small. For illustration, recall the simple example from above with

4

Figure 1.1: Blue: regularized weights {ω ∈ Ω2,2 : min (ASB(ω), W SB(ω)) ≥ δ},
orange: {ω ∈ Ω2,2 : ASB(ω) < δ}, red: {ω ∈ Ω2,2 : W SB(ω) < δ}, for
binary alphabet A = {0, 1} with Ω2,2 as in (7) and δ ∈ {0.01, 0.1, 0.2, 0.4}
(from left to right).
A = {0, 1}, m = 2, M = 1, there the condition ASB(ω), W SB(ω) > 0 is equivalent
to ω1 6= ω2 . Note that in the univariate setting M = 1 (see (Behr et al., 2017)) it
always holds that ASB(ω) ≤ W SB(ω). This is in contrast to the multivariate setting
with arbitrary M as considered here. For example,
√ when M = m and ω = Im×m
equals the identity matrix, then ASB(Im×m ) = M mini=1,...,k−1 |ai+1 − ai | > 0 =
W SB(Im×m ).
Second, we discuss conditions on F . In order to identify ω from G it is necessary
that the sources F·1 , . . . , F·m differ sufficiently much. For instance, if F·1 = . . . = F·m
then G = F·1 irrespective of ω. Here, we employ the separability condition from
Behr and Munk (2017), which provides a sufficient variability of F and guarantees
identifiability. More precisely, separability guarantees that for each i = 1, . . . , m there
exists some j = 1, . . . , n where Fji takes the second smallest alphabet value and all
other sources Fki , k 6= j, take the smallest alphabet value. As the alphabet is of the
form (4), this is equivalent to
ω1· , . . . , ωm· ∈ {G1· , . . . , Gn· }

⇔

e1 , . . . , em ∈ {F1· , . . . , Fn· },

(8)

where ei denotes the i-th unit vector. Note that for the simple example from above
with A = {0, 1}, m = 2, M = 1 this is analog to F·1 6= F·2 .
In summary, we denote a pair (ω, F ) in MABS as δ-separable if (8) holds for the
sources F and min(ASB(ω), W SB(ω)) ≥ δ. If δ > 0, in Theorem 2.2, we extend (Behr
et al., 2017, Theorem 1.4) for M = 1 to arbitrary M ∈ N and show that δ-separability
provides perturbation stability for MABS via the parameter δ, yielding exact recovery
for the sources F in a neighborhood of the mixture G, that is for (ω, F ), (ω, F ) both
δ-separable
max k(F ω)j· − (F ω)j· k < c(δ)
j=1,...,n
(9)
⇒ F = F and max kωi· − ωi· k < c(δ),
i=1,...,m

where c(δ) → 0 as δ → 0.

5

Minimax estimation The main contribution of this paper (see Section 3 and 4) is
to derive estimators for the unknown quantities F and ω, which are asymptotically
minimax optimal (up to constants) for prediction (as n, M → ∞) and estimation (as
M → ∞), where we assume the number of sources m and the alphabet A to be fixed and
known. The aim is to study the influence of all quantities on recovery as n, M → ∞,
where the total sample size is nM . It is intuitive that increasing M makes inference
on F in (1) easier (we observe more mixtures of the same sources) while inference
on ω becomes more difficult (the dimension of ω and thus the number of parameters
to be estimated increases with M ). In contrast, increasing n makes inference on F
more difficult (the dimension of F and thus the number of parameters to be estimated
increases with n) while inference on ω becomes easier (we observe more samples from
the the same mixture matrix ω). This is quantified in Theorem 3.1 (lower bound) and
Theorem 3.3 (upper bound) which provide under weak identifiability conditions the
minimax rate for the prediction error as
!
0 M
M
σ 2 m σe−c σ2
kθ̂ − F ωk2
σ 2 m σe−c σ2
+ √
. inf sup EF ω
.
+ √
,
(10)
n
nM
n
n M
n M
θ̂ F ,ω
whenever ln(n)/M → 0, M/n → 0. Further, we show that the least squares estimator
(LSE) achieves this rate. Here . and & denote inequalities up to a universal constant
which does not depend on any model parameter and c = c(m, A), c0 = c0 (m, A) > 0
are positive constants. A major consequence of (10) is that when M  ln(n) the
unknown design matrix F in the linear model (1) does not play much of a role for the
prediction error. In particular, for M  ln(n) the precise form of the alphabet A does
not influence the prediction rate. The alphabet only enters in the second term, that
is, it determines the constant of the O(ln(n)) mixtures needed to remove the influence
of the unknown design F in (1). This dependence on A does not coincide for our lower
and upper bound and it seems natural that the precise constant is not a simple feature
of A. In our simulation study in Section 6 we analyze the dependence on A further.
The exact recovery result in (9) puts the basis to relate the prediction error kθ̂−F ωk
in (10) to the estimation error via the metric
√
d ((F , ω), (F, ω)) = ( M /m) 1F 6=F + max kωi − ωi k
(11)
i=1,...,m

√

(see Theorem 4.2), where the scaling factor M /m naturally arises from the dimensionality of kωi − ω̂i k. In Theorem 4.3 we show that
 
2 
d
(F
,
ω),
(
F̂
,
ω̂)
M
σ


σ 2 c1 + √ c2 e−c σ2 . inf sup EF ω 

M
M
F̂ ,ω̂ F ,ω
(12)
0
M
σ
. σ 2 c01 + √ c02 e−c σ2 ,
M
whenever ln(n)/M → 0, M/n → 0, with c1 = c1 (m, n, A), c2 = c2 (m, n, A), c01 =
c01 (m, A), c02 = c02 (m, A) > 0 positive constants and this rate is achieved by the LSE.

6

A major consequence of (12) is that also for the estimation error if M  ln(n) the
unknown F in (1) does not play much of a role.
Mostly related to our results is (Pananjady et al., 2017), who also considered model
(1), but with F being an arbitrary design matrix (without finite alphabet constraint)
which is unknown only up to a permutation matrix. They derive minimax prediction
rates for this model and show that the LSE obtains the optimal rates (up to logfactors). They also consider the case where F is unknown up to a selection matrix
(i.e. not every row of the design necessarily appears in the data Y and some rows
might be selected several times). One can rewrite the MABS model in an analog
way, to obtain a model as in (1) where the design matrix equals F = ΠA, with Π
an unknown selection matrix and A being the matrix where the rows constitute of
all different combinations of alphabet values (see Section 2). Pananjady et al. (2017)
consider general A and derive minimax prediction rates of the form


σ 2 m σ 2 (ln(n))
1
kθ̂ − ΠAωk2 ≈
+
,
(13)
inf sup EΠAω
nM
n
M
θ̂ ΠAω
where the log-term only appears in their upper bound. In our situation, where we
assume a specific finite alphabet for the design matrix, thus a specific matrix A, the
second term in the minimax rate becomes exponential in M instead of parametric.
The rate (13) is obtained in (Pananjady et al., 2017) by treating the whole matrix ΠA
as unknown. Here we can exploit a specific structure of A and thus obtain a faster
rate. Note that, just as in our setting (see (10)), Pananjady et al. (2017) obtain with
(13) that whenever ln(n)  M the unknown permutation Π does not play much of a
role for the prediction error. The major difference, however, is that under the finite
alphabet we can now provide identifiability conditions on F = ΠA and ω in (1) and
thus, in contrast to Pananjady et al. (2017), we do obtain estimators for ω and Π and
bounds for the estimation error. By regularizing the model in an appropriate way, we
obtain the minimax estimation rate for M → ∞ for F and ω up to constants and
show that it is achieved by the LSE.
Our theoretical findings show that the LSE achieves optimal rates, both for the
prediction error and for the estimation error. Unfortunately, we are not aware of
any polynomial time algorithm for exact computation of the LSE. In fact, Pananjady
et al. (2016) have shown that the LSE for their problem (see Section 5 for details)
is NP-hard to compute. We conjecture that this is the case here, too. Therefore, in
Section 5 we propose a simple iterative Lloyd’s algorithm to approximate the LSE,
where simulations suggest similar rates as for the LSE. In particular, simulations show
that the Lloyd’s algorithm outperforms standard methods to estimate the underlying
signal F ω (not F and ω separately), which do not take the finite alphabet structure
into account, like, e.g., singular value (SV) thresholding, as proposed in (Pananjady
et al., 2017).
Related work Note that model (1) with either of both, F or ω, known corresponds to
standard models in statistics: If F is known, (1) is a classical multivariate LM. If ω is
known, the finite alphabet assumption turns model (1) into a clustering problem with

7

known centers {aω : a ∈ Am }. In contrast, in MABS both, F and ω, are unknown.
Hence, MABS can be seen as a hybrid model of clustering and parameter estimation,
but has received rather few attention so far. Only some specific instances of MABS,
e.g., specific alphabets, have been considered previously, see, for example, (Talwar
et al., 1996; Pajunen, 1997; Diamantaras, 2006; Van der Veen, 1997). However, all
these works only focus on algorithms (often in the noiseless case), but do not provide
any theory for a statistical model as in (1).
Identifiability for finite alphabet sources has been considered e.g., in (Talwar et al.,
1996; Diamantaras, 2006), see (Behr and Munk, 2017) for recent results. From there
the separability condition, which is employed in this paper, is borrowed. Here, we
extend the identifiability result of (Behr and Munk, 2017) to a stable and exact recovery
result for an arbitrary mixture dimension M ∈ N (see Theorem 2.2). The univariate
case M = 1 was considered in (Behr et al., 2017), where the temporal structure
of a change-point regression setting allows recovery of sources and mixing weights.
Here, we treat mixtures with arbitrary dimension M ∈ N and the validity of our
methodology now results from observing M (possibly different) mixtures of the same
sources. In particular, we do not require any “temporal” structure as in the univariate
case. In summary, rigorous statistical modeling, methodology, and theory for the
MABS problem as in (1) was entirely lacking, to best of our knowledge, prior to this
work.
A similar problem as in (Pananjady et al., 2017), which we discussed above, was
considered in (Pananjady et al., 2016). However, they assume M = 1 (not general M ∈
N as in this paper) and unknown permutation matrices Π, which is more restrictive
than general selection matrices, as considered here. Moreover, they assumed a random
design A with Gaussian entries, in contrast to MABS where A is a fixed finite alphabet
matrix. This makes the analysis severely different, as the finite alphabet assumption
allows to separate Π and ω from their mixture ΠAω. In their setting, they give a
sharp condition on the signal to noise ratio kωk/σ and the number of observations n
under which it is possible to exactly recover the permutation Π with large probability.
They show that the LSE recovers Π with large probability whenever this is possible.
Unnikrishnan et al. (2015) study a similar model as Pananjady et al. (2016) but mainly
focus on the noiseless case. They also consider a random design for A (in contrast to
MABS). They focus on recovery of ω (not on Π) and show that whenever n > 2m with
probability one ω can be recovered from the (noiseless) observations Y . Marques et al.
(2009) consider a similar model as Pananjady et al. (2017) in the context of object
recognition, where m = 3 and M = 2. There m = 3 corresponds to the dimension
of an object, M = 2 to the dimension of a photo of this object, and the unknown
mixture matrix ω to an unknown camera perspective. They also focus on recovery of
the unknown permutation Π. Their results basically require that sufficiently many of
the n permutations are known in advance, which is a rather strong assumption.
A structural similarity to MABS appears in nonnegative matrix factorization (NMF),
where one assumes (1) with F and ω both non-negative (Lee and Seung, 1999; Donoho
and Stodden, 2004; Arora et al., 2012). Here, however, we do have the additional assumption of a finite alphabet for the sources F , which leads to a model structure more
related to a classification problem. Hence, estimation rates are expected to be in a

8

completely different regime, although, we stress that we are not aware of any minimax
results for NMF, as derived here for MABS. Further, it should be stressed that from a
computational perspective both models are very different, as NMF algorithms usually
do not allow to incorporate a finite alphabet assumption. An exception is binary matrix factorization (Li, 2005), where both, F and ω are assumed to have binary entries.
There, however, also the data matrix Y is assumed to be binary, which makes both,
theory as well as computations very different.
The separability conditions (δ-IC 1) and (IC 2), which we employ here, are closely
related to standard conditions for NMF (Donoho and Stodden, 2004; Arora et al.,
2012), from where the notation separable originates. The proofs of these conditions
are, however, not similar. Whereas in MABS they are build on combinatorics, in
NMF they are build on geometric considerations. More details on the relation between
separability in MABS, as it is employed here, and separability in NMF can be found,
for example, in (Behr et al., 2017; Behr and Munk, 2017).
Another related problem occurs in Independent Component Analysis (ICA) (see e.g.,
(Comon, 1994)), which is based on the stochastic independence of the different sources
(assumed to be random and non-Gaussian) and which assumes i.i.d. observations from
the mixture distribution. Here we do not make any distributional assumption on the
different sources. We rather treat them as deterministic and fixed, which reveals the
MABS model (1) rather as a regression problem. Also related is statistical seriation
(Flammarion et al., 2016) which, just as the MABS model, relates to a “classical”
model with additional combinatorial uncertainty. There, however, F in (1) is itself
a permutation matrix and ω is assumed to have unimodal rows, which is a rather
different setting. Interestingly, just as in this paper for MABS, (Flammarion et al.,
2016) obtain that the LSE is (almost) minimax optimal for the statistical seriation
problem.
Organization of the paper The paper is organized as follows. In Section 2 we introduce the MABS model and corresponding identifiability conditions. From this we
derive stable recovery under suitable regularization, see Theorem 2.2. In Section 3 we
derive lower bounds (Theorem 3.1) and upper bounds (Theorem 3.3) for the minimax
prediction rate. In Section 4 we give an upper bound for the classification error (Corollary 4.1) and, based on this, lower and upper bounds for the minimax estimation rate
(Theorem 4.3). A Lloyd’s algorithm for approximating the LSE is proposed in Section
5 and simulation results are given in Section 6. We conclude in Section 7. All proofs
are postponed to the supplement.

2 Model assumptions and identifiability
It is illustrative to rewrite the MABS model to highlight its combinatorial structure.
To this end, we rewrite the unknown finite alphabet design matrix F ∈ An×m as a
product of an unknown selection matrix Π and the known design matrix A with rows
consisting of all different alphabet combinations in Am . Then the MABS model (1) is

9

equivalent to
Y = ΠAω + Z,

(14)

with an unknown selection matrix
m

Π ∈ {0, 1}n×k ,

n
X

Πij = 1 ∀i = 1, . . . , n,

(15)

j=1

and known finite alphabet design matrix

0
0
0 ... 0
0
0
0
0
.
.
.
0
0

0
0
0
.
.
.
0
0


..

.
A := 
0
0
0 ... 0
0

0
0
0 ... 0
1


..

.
ak ak ak . . . ak ak


0
1

a3 



 ∈ {0, 1, a3 , . . . , ak }km ×m ,
ak 

0



ak

(16)

where the rows of A constitute all different vectors in Am (recall the alphabets as
in (4)). Further, the unknown mixing matrix ω ∈ Ωm,M is as before in (1) and we
assume i.i.d. normal noise Zij ∼ N (0, σ 2 ), i = 1, . . . , n, j = 1, . . . , M , with unknown
variance σ 2 . As discussed in the introduction, we employ δ-separability on Π and ω
to guarantee identifiability and stable recovery. That is, for ω ∈ Ωm,M and δ > 0, we
assume the identifiability condition (IC)
min (ASB(ω), W SB(ω)) ≥ δ.

(δ-IC 1)

The separability condition for F in (8) now translates to Π in (15) as follows. Π is
separable if and only if F = ΠA is separable, namely,
∃j1 , . . . , jm ∈ {1, . . . , n} such that ((ΠA)rs )r=j1 ,...,jm = Im×m ,

(IC 2)

s=1,...,m

that is, if Π selects at least once each of the unit vectors e1 , . . . , em .
Theorem 2.1 (Behr and Munk, 2017, Theorem 4.1).
Consider the MABS model (14). Then, for any δ > 0 under (δ-IC 1) and (IC 2),
(ω, Π) is identifiable:
If (ω, Π), (ω, Π) both satisfy (δ-IC 1) and (IC 2) for some δ > 0, then ΠAω = ΠAω
implies ω = ω and Π = Π.
The following theorem shows how the parameter δ in (δ-IC 1) regularizes the space
of possible mixtures ΠAω in such a way that stable recovery is obtained within perturbations of G = ΠAω.

10

Theorem 2.2 (Exact recovery).
Consider the MABS model (14). Fix δ > 0 and let 0 <  < δ/(1 + mak ). Assume
that (δ-IC 1) and (IC 2) hold for (ω, Π), (ω, Π). If
max k(ΠAω)j· − (ΠAω)j· k < ,

j=1,...,n

1. then maxi=1,...,m kωi· − ωi· k <  and
2. Π = Π.
In words, whenever two δ-separable mixtures are close enough, the corresponding
weights are also close and the corresponding sources are even the same.
Remark 2.3 (Converse exact recovery).
Note that the converse direction of Theorem 2.2 also holds up to a constant factor.
More precisely, for any  > 0, if for some (ω, Π), (ω, Π) it holds that maxi=1,...,m kωi· −
ωi· k <  and Π = Π, then maxj=1,...,n k(ΠAω)j· − (ΠAω)j· k < mak . This follows
directly from the triangle inequality.
Remark 2.4 (Separability for general mixing weights).
For arbitrary mixing matrices, not necessarily positive and summing up to one,
ω ∈ ω ∈ Rm×M : 0 < kω1· k < . . . < kωm· k , separability (IC 2) can be defined analogously, where the uncertainty about sign and scaling of the weights leads to additional
vectors which must be selected by Π to guarantee identifiability, see (Behr and Munk,
2017, Theorem 7.1). Theorem 2.1 and 2.2 can be adapted accordingly.
Analog to δ which regularizes the identifiability condition on ω, we further introduce
a second parameter which regularizes separability of Π in (IC 2). To this end, ΠA is
denoted as Λ-separable, if the i-th unit vector ei appears in ΠA at least Λ times for
each i = 1, . . . , m, that is

# j ∈ {1, . . . , n} : (ΠA)j· = ei ≥ Λ for all i = 1, . . . , m.
(Λ-IC 2)
Thus, for fixed alphabet A, number of sources m, and A as in (16), we consider the
δ,Λ
parameter space of (δ, Λ)-regularized response matrices Nn,M
=

N δ,Λ := ΠAω : ω ∈ Ωm,M , Π as in (15) with (δ-IC 1), (Λ-IC 2) .

(17)

3 Minimax rates for the prediction error
The following theorem gives a lower bound on the minimax prediction error in the
MABS model. To this end, for the fixed alphabet A, define the smallest pairwise
difference
∆Amin := min {|a − a0 | : a 6= a0 ∈ A} .

11

(18)

Theorem 3.1 (Lower bound).
Consider the MABS model (14)
parameter space N δ,Λ as in
√ Λ)-regularized
√ with (δ,
2
(17). Further, assume that σ/ 8 < δ ≤ M (∆Amin ) (90ak m)−1 and 1 ≤ Λ ≤ bn/mc.
Let min(n, M ) ≥ 2m. Then


EΠAω kθ̂ − ΠAωk2
inf
sup
θ̂ ΠAω∈N δ,Λ


≥ 0.4M

1
288m5 a2k
+
2
(m − 1)σ
(∆Amin )2 n

−1

δ2
1
+ σδe− 8σ2 .
2

The next theorem gives an upper bound on the minimax prediction error. If δ and
Λ are chosen appropriately, it almost coincides with the lower bound from Theorem
3.1 (see Corollary 3.4). This upper bound is achieved by the LSE
θ̂ ∈ argminθ̃∈N δ,Λ kY − θ̃k2 .

(19)

Remark 3.2 (LSE).
a) (Existence) The LSE θ̂ exists, that is the minimum in (19) is attained. To see
this, note that N δ,Λ =
[
{ΠAω : ω ∈ Ωm,M , min(ASB(ω), W SB(ω)) ≥ δ}
(20)
Π Λ-separable

is a finite union of closed, bounded in [0, ak ]n×M , hence, compact sets.
b) (Uniqueness) θ̂ is not always unique, but the upper bound in Theorem 3.3 holds
for any minimizer. A counterexample is the following. Let n = m = 2, M = 1,
Λ = 1, A = {0, 1, 2}, and δ = 0.01. By separability, Π̂A is restricted to the
identity matrix. Thus
argminθ̃∈N δ,Λ kY − θ̃k2 = argminω∈Ωδ2,1 kY − ωk2
= argminω∈Ωδ2,1 ω12 − (1 + Y1 + Y2 )ω1 ,
with ω1 = 1 − ω2 . Simple calculations give that Ωδ2,1 = {(ω1 , 1 − ω1 )> : ω1 ∈
[0.1, 3/10] ∪ [11/30, 0.45]}. If the observations Y1 , Y2 are such that (1 + Y1 +
Y2 )/2 = 1/3, it is easy to check that
argminω∈Ωδ2,1 ω12 − (1 + Y1 + Y2 )ω1 = {(3/10, 7/10)> , (11/30, 19/30)> }.
c) (Computation) We are not aware of an efficient implementation of the LSE and
we speculate that this is an NP-hard problem in general (see the discussion in
Section 5). In Section 5 an iterative Lloyd’s algorithm for its approximation is
proposed.

12

Theorem 3.3 (Upper bound).
Consider the MABS model (14) with (δ, Λ)-regularized parameter space N δ,Λ as in
(17) and let θ̂ be the LSE in (19). Then


sup EΠAω kθ̂ − ΠAωk2
ΠAω∈N δ
Λδ 2
M 5/2 −
≤ 4σ 2 mM + 12σn2 k m m7/2 a3k √ e 8(mM +1)(1+mak )2 σ2 .
Λδ

Let throughout the following & and . denote inequalities up to a universal constant
which does not depend on any model parameter. Further, for some constants c1 , c2 > 0
as in Corollary 3.4 assume that


144σ 2 m3 a3k 
3 7/3 m
,
2
ln(n)
+
2
ln(M
)
+
ln
a
m
k
k
c2 c21

 6 2 2
m ak σ
, 2c2 mM .
n ≥ min
(∆Amin )2

M≥

(21)

Note that (21) holds eventually, whenever ln(n)/M → 0 and M/n → 0.
Corollary 3.4 (Minimax prediction rate).
√
Assume the setting of Theorem 3.3, such that δ √= c1 M and Λ = c2 M ≥ 1 for
some constants c1 , c2 > 0. Further, assume that σ/ 8M < c1 < (∆Amin )2 (90ak m)−1 .
Then, whenever (21) holds,


c2
1 M
1
σ2 m
σc1
EΠAω
inf
sup
kθ̂ − ΠAωk2 &
+ √ e− 8 σ 2 ,
nM
n
n M
θ̂ ΠAω∈N δ,Λ


M
1
σ2 m
σ
2
√ e−c σ2 ,
EΠAω
sup
kθ̂ − ΠAωk .
+√
inf
nM
n
c2 c1 n M
θ̂ ΠAω∈N δ,Λ
with c := c2 c21 /(16m3 (1 + ak )2 ), and the LSE achieves the second inequality.
Remark 3.5 (Separability as M →√∞).
Note that the assumption δ = c1 M and Λ = c2 M for some constants c1 , c2 in
Corollary 3.4 is not very restrictive. For √
example, when ω is uniformly
distributed,
√
the ASB(ω) in (δ-IC 1) grows with rate M , that is ASB(ω)/ M ≥ c1 for some
constant c1 asymptotically with probability one as M → ∞, see Lemma S2.2. Further,
by (Behr and Munk, 2017, Theorem 5.1), when F = ΠA is chosen such that Fij ,
i = 1, . . . , n, j = 1, . . . , m, are independent and take any value in A equally likely,
for Λ as in (Λ-IC 2) it holds that Λ/M ≥ c2 for any c2 > 0 with probability one, as
n/M → ∞.
Corollary 3.4 sheds light on the specific tradeoff between n and M regarding the
prediction error of ΠAω. The dependence on M vanishes exponentially fast (as long
as ln(n)/M → 0 and M/n → 0). Hence, for sufficiently large M the prediction rate is
dominated by its first term, which is parametric in n. Thus, as long as ln(n)/M → 0

13

and M/n → 0 the unknown selection matrix Π in the linear model (14) does not play
much of a role. Put it differently, if, in a multivariate linear model (1) the dimension
M (the number of mixtures of the BSS problem (2)) is at least of order O(ln(n)), then,
under weak identifiability conditions, for the estimation accuracy of the signal F ω it is
irrelevant whether the design F is completely known or it is only known up to a finite
set of possible values for its entries. This observation was already made for the more
general setting of arbitrary design A with unknown selection matrices Π in (Pananjady
et al., 2017). However, here the additional assumption of the given finite alphabet
constraint leads to a much faster decay (recall our discussion in the introduction).
More importantly, the finite alphabet assumption provides identifiability and hence,
allows separation of ω and Π, such that the results on the prediction error can be
extend to the estimation error, as we will show in the following. This is not possible
for the more general setting as in (Pananjady et al., 2017).

4 Estimation error
From the proof of Theorem 3.3 it is easy
 to derive the following upper bound on the
maximal classification error P Π̂ 6= Π .
Corollary 4.1 (Upper bound on classification error).
Assume the setting of Corollary 3.4 and let θ̂ be the LSE in (19). Then, whenever
(21) holds
2

sup
ΠAω∈N δ,Λ



PΠAω Π̂ 6= Π . √

c2 c1
σ
−
√ e 16m3 (1+ak )2
c2 c1 n M

M
σ2

.

(22)

In order to derive lower bounds for the maximal estimation error, one can combine
Corollary 3.4 with Theorem 2.2. To this end, recall the metric d in (11) with F replaced
by Π, which combines the classification Π̂ 6= Π and estimation kωi − ω̂i k error. The
metrics d ((Π, ω), (Π, ω)) and kΠAω − ΠAωk are locally equivalent on N δ,1 as the
following theorem shows.
Theorem 4.2.
Let ΠAω, Π0 Aω 0 ∈ N δ,1 as in (17), then
1.

√
d ((Π, ω), (Π0 , ω 0 )) ≥ kΠAω − Π0 Aω 0 k/( nmak ),

2. if kΠAω − Π0 Aω 0 k ≤ δ/(1 + mak ), then
d ((Π, ω), (Π0 , ω 0 )) ≤ kΠAω − Π0 Aω 0 k.
The following theorem shows that the LSE is not only asymptotically minimax rate
optimal for the prediction error as in Corollary 3.4, but also as M → ∞ minimax rate
optimal for the estimation error in terms of the metric d in (11) .

14

Theorem 4.3 (Minimax estimation rate).
Assume the setting of Corollary 3.4. Then, whenever (21) holds
 
2 
c2
√
1 M
c1
1
inf sup EΠAω d (Π, ω), (Π̂, ω̂)
+σ M
e− 8 σ 2 ,
& σ2 M
2
2
2
nmak
nm ak
Π̂,ω̂ Π,ω
 
2 
3 2
√
m2 a2k −c M2
m ak
inf sup EΠAω d (Π, ω), (Π̂, ω̂)
M
+
σ
e σ ,
. σ2 M
√
c21
c2 c31
Π̂,ω̂ Π,ω
with c := c2 c21 /(16m3 (1 + ak )2 ), and the LSE achieves the second inequality.
Again, Theorem 4.3 shows that when M is sufficiently large, increasing M further
does not influence the estimation rate in terms of d(·, ·)2 /M . Moreover, for M large
enough, the minimax estimation rate of d(·, ·)2 /M does not increase with n, although
the dimension of Π and Π̂, respectively, increase with n. Thus, if ln(n)  M the
unknown selection Π in the linear model (14) does not play much of a role for the
estimation rate, as well.

5 A clustering based algorithm
In Section 3 and 4 we have seen that the LSE attains (almost) optimal rates for both,
the maximal prediction error and the maximal estimation error (in terms of the metric
d in (11)). Unfortunately, we are not aware of any polynomial time algorithm for
its exact computation. To this end, note that in (20) for Π 6= Π0 , both Λ-separable,
Theorem 2.2 implies that the corresponding two sets in the union in (20) are disjoint.
Thus, computation of the LSE amounts to minimization over exponentially many (in
n) disjoint, compact sets (see also Figure 1.1). Pananjady et al. (2016) have shown
that exact computation of the LSE is NP-hard in general, for the MABS model (14)
with M = 1, but for arbitrary design A (not the specific form in (16)) and restricted to
permutation matrices Π (not the bigger class of selection matrices). Although, their
results do not directly apply to our setting, it is near at hand that exact computation
of the LSE for the MABS model (14) is also not feasible.
Therefore, in the following a Lloyd’s algorithm for its approximation is proposed.
To ensure minimax optimality in N δ,Λ the LSE in (19) depends on the regularization
parameters (δ, Λ). For computational purposes consider now the unrestricted case δ =
Λ = 0. Then, on the one hand, given the selection Π computation of the unrestricted
LSE corresponds to a convex optimization problem1 , which can be solved easily with
standard techniques, see e.g. (Van den Meersche et al., 2009). On the other hand, given
the mixture matrix ω computation of the LSE corresponds to a simple LS clustering
with known centers Aω ∈ RK×M . Therefore, an iterative Lloyd’s algorithm (see e.g.,
(Lu and Zhou, 2016)) which starts with a random ω̂ and successively updates Π̂ and
ω̂ can be employed, see Figure 5.1 for details. Note that this algorithm may as well be
1 Note

that the condition 0 < kω1· k < . . . < kωm· k in (3) can be neglected in this optimization step,
as for any permutation matrix P it holds that kY − ΠAωk = kY − ΠAP −1 P ωk and hence, one
can reorder the columns of ω̂ and the rows of Π̂A subsequently.

15

Input: Y , A, t, m
ω̂ 1 ← uniform random choice
Π̂j· A ← argmina∈Am kYj − aω̂ 1 k
for j = 1, . . . , n
2
2
ω̂ ← argminω kY − Π̂Aωk
1
2
> t do
while max1≤i≤m, 1≤j≤n ω̂ij
− ω̂ij
1
2
ω̂ ← ω̂
Π̂j· A ← argmina∈Am kYj − aω̂ 2 k
for j = 1, . . . , n
ω̂ 2 ← argminω kY − Π̂Aωk2
end while
Return: ω̂ 1 and Π̂.
Figure 5.1: Lloyd’s algorithm for approximation of the LSE in the MABS model.
applied to the more general setting of an arbitrary design A with unknown selection
matrix Π as in (Pananjady et al., 2017).
Clearly, in each step of the Lloyd’s algorithm kY − Π̂Aω̂k2 decreases and thus,
converges to a local minimum of (ω, Π) 7→ kY − ΠAωk2 . In order to reach the global
minimum, we propose to compute K realizations of the algorithm in Figure 5.1 (for
several random initial values for ω̂) and chose the one which minimizes (ω, Π) 7→
kY − ΠAωk2 . In practice, t = 0.001 (see Figure 5.1) and K = 3 worked well, usually
terminating within 10 iterations. For instance, computation in R for n = M = 1, 000,
m = 3, A = {0, 1} takes about 5 seconds on a desktop computer with intel core i7
processor.

6 Numerical results
In the following, the Lloyd’s algorithm from Figure 5.1 is explored in a simulation
study. In particular, we want to compare its performance with the theoretical findings
for the LSE from Section 3, which cannot be computed efficiently (see Section 5 and
Remark 3.2). Corollary 3.4 yields that the LSE achieves optimal rates for the maximal
prediction error. The maximal prediction error cannot be simulated efficiently, as
the maximum may be attained at any (ω, Π). Instead, we simulate the Bayes risk
with uniform priors for ω√and Π, respectively. As discussed in Remark 3.5, this
ensures that ASB(ω) ≥ c1 M and Π is c2 M -separable, for some constant c1 , c2 > 0
asymptotically almost surely2 . We simulate the estimator θ̂ = Π̂Aω̂ (see Section 5)
of Lloyd’s algorithm from Figure 5.1 for m = 2, 3, 4 number of sources, alphabets
A = {0, 1}, {0, 1, 2, 3}, n = 500, 1000 number of observations, for M ∈ {1, . . . , 1000},
and standard deviation σ = 0.2, 0.5, 1. Simulation runs were always 100, 000.

√
uniform prior for ω does not guarantee W SB(ω) ≥ c1 M asymptotically almost surely. However, this is more of a technical issue, as one can always consider estimation up to permutation matrices P (recall that F P P −1 ω = F ω), that is one considers the equivalence class
(F , ω) ∼ (F P, P −1 ω) for any permutation matrix P .

2A

16



Figure 6.1: Solid lines: Left: Normalized MSE E kθ̂ − ΠAωk2 /M for the estimator
θ̂ = Π̂Aω̂ of Figure 5.1 for m = 3 sources, a binary alphabet A = {0, 1}, and
n = 500 observations, for M ∈ {1, . . . , 1000} and σ = 0.2, 0.5, 1 (lightgray,
darkgray, and black line).
 Top right: MSE on a logarithmic scale. Bottom
right: E kθ̂ − ΠAωk2 /(M σ 2 ). Dashed lines: results for SV thresholding
as proposed in (Pananjady et al., 2017).
Dependence on σ Simulations of the mean squared error (MSE) for m = 3, A =
{0, 1}, n = 500, M ∈ {1, . . . , 1000}, and σ = 0.2, 0.5, 1, are shown in Figure 6.1 (solid
lines). For larger variances (σ = 0.5, 1) the top left plot of Figure 6.1 (solid lines)
shows a peak at M = 3 and for M > 3 an exponential decay to some limiting value.
For small variance σ = 0.2 the exponential decay starts already at M = 1. The
peak can be explained as follows: When the variance is large, i.e., the signal to noise
ratio is small, the observations Y in (14) are likely to lie outside of the parameter
space N ⊂ [0, 1]n×M . Thus, the Lloyd’s algorithm prefers parameters in N which
are close to the boundary of [0, 1]n×M . This corresponds to mixing matrices ω with
columns equal to unit vectors, which in general are not close to the true underlying
mixing matrix and hence, lead to a larger MSE. However, if M < m such ω have
zero alphabet separation boundary ASB(ω) = 0 and in particular, small number of
possible mixture values. More precisely, for M < m when the columns of ω are unit
vectors and ω̃ is a different mixing matrix with ASB(ω̃) > 0 then
# {eω : e ∈ Am } ≤ k M < k m = # {eω̃ : e ∈ Am } .

17

Consequently, when M < m, mixing matrices ω with columns being unit vectors
are not (wrongly) selected by the Lloyd’s algorithm and hence, the estimate fits the
true underlying parameters better. On the other hand, as the (maximal) number of
mixture values k m does not grow with M , when M  m choosing the columns of ω
as unit vectors only explores certain small (relatively decreasing with M ) areas of the
boundary of the parameter space N and thus does not lead to a better fit of the data,
in general. Hence, we find that M ≈ m corresponds to the most difficult scenario for
the Lloyd’s algorithm to recover the signal. However, note that this phenomena does
not admit the following interpretation. One might think that for M ≈ m estimating
each column of ω separately from the corresponding column of the data matrix Y leads
to a smaller MSE than estimating the whole matrix ω at once from the whole data
Y . This is not true in general, which can be seen by the following example. When
M = 1 estimation of ω ∈ Rm×1 reduces to estimation of its entries ω11 , . . . , ωm1 ,
as their ordering is determined via the relation ω11 < . . . < ωm1 . However, when
M = 2 estimating its entries ωij does not determine ω uniquely from the relation
kω1· k < . . . < kωm· k. In particular, combining estimates of single columns of ω is not
feasible. The top right plot of Figure 6.1 (solid line) shows the
logarithmic
 MSE on a 
scale (where we subtracted the limiting value limM →∞ E kθ̂ − ΠAωk2 /M ). Its
linearity suggests an exponential decay confirming the bounds in Corollary 3.4. As in
Corollary 3.4 the slope in Figure 6.1 (top right, solid line) decreases with σ 2 and the
intercept does not depend on σ 2 . The bottom left plot of Figure 6.1 (solid line) shows
the limiting value of the MSE. In the bottom right plot of Figure 6.1 (solid line) one
observes that, just as in Corollary 3.4, it scales with σ 2 .

Dependence on n Simulations for m = 3, A = {0, 1}, M ∈ {1, . . . , 1000}, and
n = 500, 1000 (and σ = 0.2, 0.5) are shown in Figure 6.2. The top right plot of Figure
6.2 shows an exponential decay. As in Corollary 3.4 the slope does not depend on
n and the intercept increases with n. The bottom left plot of Figure 6.2 shows the
limiting value of the MSE, just as in Corollary 3.4 it does not depend on n.
Dependence on m Simulations for A = {0, 1}, n = 500, M ∈ {1, . . . , 1000}, and
m = 2, 3, 4 (and σ = 0.2, 0.5) are shown in Figure 6.3. As in the upper bound of
Corollary 3.4 the slope in Figure 6.3 top right decreases with m and the intercept does
not depend on m. The bottom left plot of Figure 6.3 shows the limiting value of the
MSE, just as in Corollary 3.4 it increases with m (scaling with m − 1).
Dependence on A Simulations for m = 3, n = 500, M ∈ {1, . . . , 1000}, and A =
{0, 1}, {0, 1, 2, 3} (and σ = 0.2, 0.5) are shown in Figure 6.4. The slope in Figure 6.4
top right does not depend on the alphabet. This suggests that in Corollary 3.4 the
additional (1 + mak )2 term in the upper bound is not necessary. The intercept in
Figure 6.4 top right increases with ak as in Corollary 3.4. The bottom left plot of
Figure 6.4 shows the limiting value of the MSE. Just as in Corollary 3.4 it does not
depend on the alphabet A.

18

Figure 6.2: As solid lines in Figure 6.1 with m = 3, A = {0, 1}, n = 1000, 500 (solied
and dashed line), and σ = 0.2, 0.5 (gray and black line).

Figure 6.3: As Figure 6.2 with A = {0, 1}, n = 500, m = 2, 3, 4 (dotted, dashed, and
solid line), and σ = 0.2, 0.5 (gray and black line).

19

Figure 6.4: Same as Figure 6.2 with m = 3, n = 500, A = {0, 1}, {0, 1, 2, 3} (dashed
and solid line), and σ = 0.2, 0.5 (gray and black line).
Estimation error Figure 6.5 shows simulation results for the estimation error for m =
3, A = {0, 1}, n = 500, M ∈ {1, . . . , 1000} and σ = 0.2, 0.5, 1. Again,
 an exponential



decay to some limiting value as M → ∞ for E d((Π̂, ω̂), (Π, ω))2 /M , P Π 6= Π̂ ,



E maxi=1,...,m kωi· − ω̂i· k2 /M , and E kΠA − Π̂Ak22 is observed. Whereas, as M →



∞, P Π 6= Π̂ → 0, the limiting value of E maxi=1,...,m kωi· − ω̂i· k2 /M and, thus,


also of E d((Π̂, ω̂), (Π, ω))2 /M , scales with σ 2 , which is in accordance with Theorem
4.3.
Comparison with singular value thresholding We also simulated the MSE of the
procedure proposed in (Pananjady et al., 2017), which corresponds to a SV thresholding and which does not take into account the specific finite alphabet structure. As
the theory suggests, this leads to a worse performance. In particular, the log-plot in
the upper right plot from Figure 6.1 shows that with a simple SV thresholding the
MSE does not decrease exponentially, in contrast to the Lloyd’s algorithm, which takes
advantage of the finite alphabet (recall the discussion in the introduction and equation
(13)). Moreover, we stress that SV thresholding does not provide estimates for Π and
ω, only for the total signal ΠAω, in contrast to the Lloyd’s algorithm from Figure 5.1.

20





Figure 6.5: From top to bottom: E d((Π̂, ω̂), (Π, ω))2 /M ; P Π 6= Π̂ (left) and



E kΠA − Π̂Ak22 (right); E maxi=1,...,m kωi − ω̂i k2 /M for θ̂ = Π̂Aω̂ as
in Figure 5.1 for m = 3, A = {0, 1}, n = 500, and σ = 0.2, 0.5, 1 (lightgray,
darkgray, and black line).

21

7 Conclusion and discussion
In this paper we introduced the Multivariate finite Alphabet Blind Separation (MABS)
model with Gaussian noise, where we imposed regularity conditions to ensure identifiability of the model parameters. Depending on these quantities, we derived lower and
upper bounds (attained by the LSE) of the maximal prediction error which coincide up
to constants. In particular, our results reveal that, due to the finite alphabet structure,
minimax rates are significantly improved, in the sense that the classification error vanishes exponentially as the number of mixtures grows, instead of parametric estimation
rates of convergence which are optimal without the finite alphabet assumption (recall
the discussion in Section 1 regarding the results of (Pananjady et al., 2017)). Most importantly, we could derive bounds for the estimation error from those for the prediction
error. This only becomes feasible due to the finite alphabet structure, which provides
identifiability and stable recovery under reasonable conditions. Again, we showed that
the LSE attains the optimal rates for the estimation error. In particular, our results
demonstrate that the unknown design does not influence the minimax rates when the
number of mixtures M is at least of order ln(n), where n is the number of observations.
This is in strict contrast to a computational view on the MABS model. Whereas for
known design computation of the LSE amounts to a convex optimization problem, for
unknown finite alphabet design as in (1) it amounts to minimization over a disjoint
union of exponentially many (in n) sets. Therefore, we propose a simple Lloyd’s algorithm for its approximation. Simulations indicate similar convergence properties as
in the theoretical results for the LSE. In particular, our simulations confirm a significantly better convergence rate as methods which do not take into account the finite
alphabet structure, such as SV thresholding.
In summary, a major consequence of this paper is that finite alphabet structures can
significantly improve prediction accuracy and provide identifiability, even enabling us
to estimate Π and ω, separately. In a broader context, finite alphabet structures may
be considered as a new type of sparsity – promising to be explored further.
For future research it would be interesting to study the estimation errors of the
mixing matrix ω and the selection matrix Π more detailed. So far, we have focused on
the prediction error of their product and our corresponding results for the estimation
error are stated in terms of a metric which depends on both, the estimation error of
ω and of Π. We speculate that it might be possible to estimate one of them while
estimation of the other is not feasible. An example where such a situation occurs is
multi-reference alignment (Bandeira et al., 2017). There, one observes Yi = Ri θ+σ 2 Zi ,
where θ ∈ Rd is an unknown parameter, Ri are unknown cyclic shifts, Zi is i.i.d.
standard Gaussian noise, and σ 2 is a known variance. In this setting it is shown in
(Bandeira et al., 2017) that the parameter θ can be estimated consistently from the
data Y even when the noise level is so high that it is not possible to estimate the
shifts Ri . Analogously, it might be possible to estimate ω from Y = ΠAω + Z in (14)
without estimating Π or vice versa.

22

Acknowledgements
Merle Behr acknowledges support of RTG 2088 and CRC 803 Z2 (German Research
Foundation). Axel Munk acknowledges support of CRC 803 Z2 (German Research
Foundation). Helpful comments of Martin Wainwright, Boaz Nadler, and Burkhard
Blobel are gratefully acknowledged.

References
Arora, S., Ge, R., Kannan, R., and Moitra, A. (2012). Computing a nonnegative matrix
factorization–provably. In Proceedings of the Forty-Fourth Annual ACM Symposium
on Theory of Computing, pages 145–162. New York.
Bandeira, A., Rigollet, P., and Weed, J. (2017). Optimal rates of estimation for multireference alignment. arXiv preprint arXiv:1702.08546.
Behr, M., Holmes, C., and Munk, A. (2017). Multiscale blind source separation. arXiv
preprint arXiv:1608.07173. To appear in The Annals of Statistics.
Behr, M. and Munk, A. (2017). Identifiability for blind source separation of multiple finite alphabet linear mixtures. IEEE Transactions on Information Theory,
63(9):5506–5517.
Beroukhim, R., Mermel, C. H., Porter, D., Wei, G., Raychaudhuri, S., Donovan,
J., et al. (2010). The landscape of somatic copy-number alteration across human
cancers. Nature, 463(7283):899–905.
Carter, S. L., Cibulskis, K., Helman, E., McKenna, A., Shen, H., Zack, T., et al.
(2012). Absolute quantification of somatic DNA alterations in human cancer. Nature
Biotechnology, 30(5):413–421.
Comon, P. (1994). Independent component analysis, a new concept? Signal Processing,
36:287–314.
Diamantaras, K. I. (2006). A clustering approach for the blind separation of multiple
finite alphabet sequences from a single linear mixture. Signal Processing, 86(4):877–
891.
Donoho, D., Liu, R. C., and MacGibbon, B. (1990). Minimax risk over hyperrectangles
and implications. The Annals of Statistics, 18(3):1416–1437.
Donoho, D. and Stodden, V. (2004). When does non-negative matrix factorization give
a correct decomposition into parts? In Advances in Neural Information Processing
Systems (NIPS) 16. MIT Press, Cambridge.
Flammarion, N., Mao, C., and Rigollet, P. (2016). Optimal rates of statistical seriation.
arXiv preprint arXiv:1607.02435.

23

Greaves, M. and Maley, C. C. (2012).
481(7381):306–313.

Clonal evolution in cancer.

Nature,

Ha, G., Roth, A., Khattra, J., Ho, J., Yap, D., Prentice, L. M., et al. (2014). TITAN:
Inference of copy number architectures in clonal cell populations from tumor wholegenome sequence data. Genome Research, 24(11):1881–1893.
Jiang, Y., Qiu, Y., Minn, A. J., and Zhang, N. R. (2016). Assessing intratumor
heterogeneity and tracking longitudinal and spatial clonal evolutionary history by
next-generation sequencing. Proceedings of the National Academy of Sciences,
113(37):E5528–E5537.
Lee, D. and Seung, S. (1999). Learning the parts of objects by non-negative matrix
factorization. Nature, 401:788–791.
Li, T. (2005). A general model for clustering binary data. In Proceedings of the
Eleventh ACM SIGKDD International Conference on Knowledge Discovery in Data
Mining, pages 188–197. Chicago.
Liu, B., Morrison, C. D., Johnson, C. S., Trump, D. L., Qin, M., Conroy, J. C., et al.
(2013). Computational methods for detecting copy number variations in cancer
genome using next generation sequencing: Principles and challenges. Oncotarget,
4(11):1868.
Love, D. J., Heath, R. W., Lau, V. K., Gesbert, D., Rao, B. D., and Andrews, M.
(2008). An overview of limited feedback in wireless communication systems. IEEE
Journal on Selected Areas in Communications, 26(8):1341–1365.
Lu, Y. and Zhou, H. H. (2016). Statistical and computational guarantees of Lloyd’s
algorithm and its variants. arXiv preprint arXiv:1612.02099.
Marques, M., Stošić, M., and Costeira, J. (2009). Subspace matching: Unique solution to point matching with geometric constraints. In IEEE Twelfth International
Conference on Computer Vision, pages 1288–1294. Kyoto.
Pajunen, P. (1997). Blind separation of binary sources with less sensors than sources.
In International Conference on Neural Networks, pages 1994–1997. Houston.
Pananjady, A., Wainwright, M. J., and Courtade, T. A. (2016). Linear regression with
an unknown permutation: Statistical and computational limits. In 54th Annual
Allerton Conference on Communication, Control, and Computing. IEEE, Monticello, IL.
Pananjady, A., Wainwright, M. J., and Courtade, T. A. (2017). Denoising linear
models with permuted data. arXiv preprint arXiv:1704.07461.
Proakis, J. G. (2007). Digital Communications. McGraw-Hill series in electrical and
computer engineering. McGraw-Hill, Boston.

24

Sampath, H., Stoica, P., and Paulraj, A. (2001). Generalized linear precoder and
decoder design for MIMO channels using the weighted MMSE criterion. IEEE
Transactions on Communications, 49(12):2198–2206.
Shah, S. P., Roth, A., Goya, R., Oloumi, A., Ha, G., Zhao, Y., et al. (2012). The
clonal and mutational evolution spectrum of primary triple-negative breast cancers.
Nature, 486(7403):395–399.
Talwar, S., Viberg, M., and Paulraj, A. (1996). Blind separation of synchronous cochannel digital signals using an antenna array. I. Algorithms. IEEE Transactions
on Signal Processing, 44(5):1184–1197.
Unnikrishnan, J., Haghighatshoar, S., and Vetterli, M. (2015). Unlabeled sensing with
random linear measurements. arXiv preprint arXiv:1512.00115.
Van den Meersche, K., Soetaert, K., and Van Oevelen, D. (2009). Xsample(): An
R function for sampling linear inverse problems. Journal of Statistical Software,
30:1–15.
Van der Veen, A.-J. (1997). Analytical method for blind binary signal separation.
IEEE Transactions on Signal Processing, 45(4):1078–1082.
Verdu, S. (1998). Multiuser Detection. Cambridge University Press, Cambridge.
Yau, C., Papaspiliopoulos, O., Roberts, G. O., and Holmes, C. (2011). Bayesian nonparametric hidden Markov models with applications in genomics. Journal of the
Royal Statistical Society: Series B (Statistical Methodology), 73(1):37–57.
Zhang, Y. and Kassam, S. A. (2001). Blind separation and equalization using fractional
sampling of digital communications signals. Signal Processing, 81(12):2591–2608.

25

Supplement to
Minimax estimation in linear models with
unknown finite alphabet design
Merle Behr∗ and Axel Munk∗,†
Institute for Mathematical Stochastics∗ ,
University of Göttingen,
and
Max Planck Institute for Biophysical Chemistry† ,
Göttingen, Germany
Email: {behr,munk}@math.uni-goettingen.de

S1 Proof of the exact recovery result Theorem 2.2
Proof of Theorem 2.2. The separability condition implies that there exists ei , ẽi ∈ Am
for i = 1, . . . , m such that
kωi· − ẽi ωk <  and kωi· − ei ωk

< .

(23)

We start with proving the first assertion by induction for i = 1, . . . , m. If either e1 or
ẽ1 equals the unit vector (1, 0, . . . , 0) ∈ Rm , (23) yields
kω1· − ω1· k < .

(24)

If e1 or ẽ1 equals the zero vector (0, . . . , 0) ∈ Rm , then ASB(ω), ASB(ω) ≥ δ and
(23) contradict. So assume that e1 and ẽ1 both neither equal the first unit vector nor
the zero vector and, in particular,
e11 ≤ 1 ⇒

m
X

e1i ≥ 1

(25)

i=2

and analog for ẽ1 . W.l.o.g. assume that kω1· k ≥ kω1· k. Then
!2
M
m
X
X

e1i ωij
≥ min a22 kω1· k2 , kω2· k2 ,
ke1 ωk2 =
j=1

(26)

i=1

where the inequality follows from separating into the following cases. If e11 ≥ a2 > 1,
then
!2
m
X
1
2
2
ei ωij
≥ (e11 )2 ω1j
≥ a22 ω1j
.
i=1

If

e11

≤ 1, then by (25) ∃r > 1 such that e1r ≥ 1, and
!2
m
X
1
2
2
ei ωij
≥ (e1r )2 ωrj
≥ ωrj
.
i=1

1

(27)

In particular, (26) gives
ke1 ωk − kω1· k ≥ min (a2 kω1· k − kω1· k, kω2· k − kω1· k)
2δ
= min (k(a2 − 1, 0, . . . , 0)ωk, kω2· k − kω1· k) ≥
>
(1 + mak )

(28)

and (23) gives
ke1 ωk − kω1· k ≤ ke1 ωk − kω1· k ≤ ke1 ω − ω1· k < .

(29)

(28) and (29) contradict, which shows (24).
Now assume that
kωi· − ωi· k <  for i = 1, . . . , r − 1

(30)

and w.l.o.g. assume that
kωr· k ≥ kωr· k.
First, assume that

Pm

r
i=r+1 ei

(31)

≥ 1. Then it follows from (23) that

kωr· k =kωr· − er ω + er ωk ≥ ker ωk − kωr· − er ωk ≥ ker ωk − 
≥kωr+1· k −  ≥ kωr· k + 2δ/(1 + mak ) −  > kωr· k,

(32)

where for the third inequality we used an analog argument as in (27). (32) contradicts
(31). Thus, it follows that
err+1 = . . . = erm = 0.
(33)
Further, if err = 0, then
kωr· − er ωk ≤ kωr· − er ωk + ker ω − er ωk
≤ kωr· − er ωk + (r − 1)ak  ≤ (1 + (r − 1)ak ) < δ,

(34)

where the second inequality follows from (30) and third inequality from (23). ( 34)
and ASB(ω) ≥ δ contradict. Thus, it follows that
err ≥ 1.

(35)

Note that (33), (35) and (23) imply that
r

 > kωr· − e ωk = kωr· −

ωr· −

(err

− 1)ωr· +

r−1
X

!!
eri ωi·

k

i=1

≥ kωr· − (ωr· − x) k,
with xj =


 (er − 1)ωr· + Pr−1 er ωi·
r

ω

rj

i=1

i

if



j

(err − 1)ωr· +

otherwise.

2

Pr−1
i=1

eri ωi·


j

≤ ωrj ,

(36)

As x, ωr , and ωr − x have non-negative entries, it also follows that
kωr k ≥ kωr − xk.
(37)
Pm
Pr−1
Now, assume that i=1 eri ≥ 1. If i=r+1 ẽri ≥ 1 then it follows from (23), (36), and
(37) that
kωr+1 k − kωr k ≤ kẽr ωk − kωr − xk
≤ k (ẽr ω − ωr ) + (ωr − (ωr − x)) k ≤ 2 < 2/(1 + mak )δ,
which contradicts W SB(ω) ≥ δ and hence, it follows that ẽrr+1 = . . . = ẽrm = 0.
Pr−1
Further, if ẽrr = 0, then ẽrr ω = i=1 ẽri ωi· and
kωr − ẽr ωk ≤ kωr − ẽr ωk + kẽr ω − ẽr ωk ≤  + (r − 1)ak  < δ,

(38)

where for the second inequality we used (23) and (30). (38) and ASB(ω) ≥ δ contradict. Thus it follows that ẽrr ≥ 1. However, this implies that
r

kωr − ẽ ωk ≥ k

r−1
X
i=1

≥k

r−1
X

ẽri ωi

+

r−1
X

eri ωi + (err − 1)ωr + (ẽrr − 1)ωr k − kωr − er ωk

i=1

eri ωi k −  ≥ δ −  > ,

i=1

which contradicts (23). Hence, it follows that
er1 = . . . = err−1 = err+1 = . . . = erm = 0 and err ≥ 1.

(39)

Thus, by (31) and (23) it follows that
kωr − ωr k ≤ kerr ωr − ωr k = ker ω − ωr k ≤ ,

(40)

where for the first inequality we used the fact that for two vectors a, b with kak ≥ kbk
and a constant c ≥ 1 it follows that kca − bk ≥ ka − bk. Thus, the first assertion follows
by induction.
To show the second assertion, assume the contrary, i.e., that there exists e 6= e0 ∈ Am
such that
 > keω − e0 ωk ≥ keω − e0 ωk − ke0 ω − e0 ωk.
(41)
As ASB(ω) ≥ δ, it follows that keω − e0 ωk ≥ δ and by the first assertion of the
theorem, it follows that ke0 ω − e0 ωk ≤ mak . Therefore (41) implies  > δ − mak ,
which is a contradiction.

S2 Additional lemmas on the ASB
For the following considerations we define the space of δ-separable mixing weights as
Ωδm,M := {ω ∈ Ωm,M : min (ASB(ω), W SB(ω)) ≥ δ} .

3

(42)

Lemma S2.1.
√
If Ωδm,M in (42) is non-empty for some m, M ∈ N, then δ ≤ √M (1+mak ) .
2m(m+1)

Proof. If Ωδm,M in (42) is non-empty, then there exists an ω ∈ Ωm,M with
v
uM
uX
2
δ ≤ ASB(ω) ≤ t
ω1j
j=1

and
δ ≤ W SB(ω) ≤ (1 + mak )/2 (kωi· k − kωi−1· k)
p
≤ (1 + mak )/2 kωi· k2 − kωi−1· k2
for all i = 1, . . . , m, with ωmj = 1 − ω1j − . . . − ωm−1j . In particular, there exists
ω ∈ Rm×M
with
+
M
M
X
X

4δ 2
2
2
2
≤
min
ω
,
min
ωij
− ωi−1j
,
1j
2
i=2,...,m−1
(1 + mak )
j=1
j=1
M
X

2
(1 − ω1j − ... − ωm−1j )2 − ωm−1j



.

j=1

Moreover,
(1 − ω1j − . . . − ωm−1j )2
= 1 + (ω1j + . . . + ωm−1j )2 − 2ω1j − . . . − 2ωm−1j
=1−

m−1
X

ωij (2 −

m−1
X

2
2
− . . . − ωm−1j
.
ωsj ) ≤ 1 − ω1j

s=1

i=1

|

{z

≥1≥ωij

}

And thus,
δ2 4
(1 + mak )2
M
X
2
≤ min
ω1j
,
j=1

min

i=2,...,m−1

M
X

M
 X

2
2
2
2
ωij
− ωi−1j
,
1 − ω1j
− . . . − 2ωm−1j

j=1

j=1

≤ max
min (x1 , x2 , . . . , xm−1 , (M − mx1 − . . . − 2xm−1 )) =
m−1
x∈R

where for the second inequality we used x1 :=
PM 2
j=1 ωi−1j for i = 2, . . . m.

4

PM

j=1

2M
,
m(m + 1)

2
ω1j
and xi :=

PM

j=1

2
ωij
−

For the following considerations, we extend the definition in (18) for the fixed alphabet A and number of sources m by
∆Am := {e1 − e2 : e1 6= e2 ∈ Am } ,
2

m

∆ A

m

:= {e1 − e2 : e1 6= e2 ∈ ∆A } ,

∆Amin := min 1 |x| ,
x∈∆A

2

∆ Amin := min
|x| .
2 1

(43)

x∈∆ A

Define the constants
√

2 ∆2 Amin
c = c(m, A) := √
,
2m
3k m2 (m − 1)

√
2(1 + mak )
.
C = C(m, A) := p
m(m − 1)

(44)

Lemma S2.2.
If ω is uniformly distributed on Ωm,M in (3), then for c, C as in (44) it holds almost
surely that
ASB(ω)
ASB(ω)
c < lim inf √
≤ lim sup √
< C.
M →∞
M
M
M →∞
√
Proof of Lemma S2.2. It follows from
p Lemma S2.1 that ASB(ω)/ M is bounded
(surely) from above by (1 + mak )/ 2(m(m − 1)), which shows the inequality on the
right hand side. Further, if M = 1 and ω is drawn uniformly,
√ then it can be shown
that P(ASB(ω) > δ) ≥ 1 − d δ with d = k 2m m2 (m − 1)/( 2∆2 Amin ). For arbitrary
M ∈ N, if ω is drawn uniformly its ASB
PM is bounded by the corresponding ASB’s of the
single components, i.e., ASB(ω)2 ≥ j=1 ASB(ω·j )2 , where ASB(ω·j ), j = 1, . . . , M ,
are independent and identically distributed with
Z ∞

√
1
E ASB(ω·j )2 ≥
(1 − d x)+ dx = 2 .
3d
0
Hence, for c <

lim inf
M →∞

√1
3d

it follows from the strong law of large numbers that almost surely

M

1 X
ASB(ω)2
≥ lim inf
ASB(ω·j )2 = E ASB(ω·j )2 > c2 ,
M →∞ M
M
j=1

which shows the inequality on the left hand side.
Lemma S2.3.
min
Ω0.2∆A
as in (42) is non-empty for any m ∈ N, with ∆Amin as in (43).
m,m
Proof. For 1/2 > δ > 0 define


m−1
0

0
m
−2

2δ

δ
ω := Im×m −

1 + mak 

0
0
−(m − 1) −(m − 2)

5

...
...
..
.

0
0

...
...

1
−1


0
0


,

0
0

(45)

where Im×m denotes the m × m identity matrix. As
i = m − 1, . . . , 1 that
1−

2δ
1+mak (m

− 1) < 1, it holds for

2δ
2δ
2δ
(i − 1) − 1 −
i ≥
1 + mak
1 + mak
1 + mak

and thus W SB(ω δ ) ≥ δ. Consequently, if ASB(ω δ ) ≥ δ it follows that Ωδm,m is
non-empty. We have that
ASB(ω δ ) ≥ ASB(Im×m )
p
2δ
max m ((m − 1)(e1 − em ))2 + . . . + (em−1 − em )2
−
1 + mak e6=0∈∆A
p
2δ
= ∆Amin −
max m ((m − 1)(e1 − em ))2 + . . . + (em−1 − em )2
1 + mak e6=0∈∆A
v
r
um−1
uX
2δ
δ4a
m(m − 1)
k
≥ ∆Amin −
2ak t
,
(m − i) = ∆Amin −
1 + mak
1
+
ma
2
k
i=1
which implies that if
p
√
δ 8ak m(m − 1)
δ ≤ ∆Amin −
1 + mak

(46)

then Ωδm,m is non-empty. As
√
1+

!−1
p
8ak m(m + 1)
≥
1 + mak

r
1+

8(m + 1)
m

!−1
≥ 0.2

(46) holds for δ = 0.2∆Amin .
Lemma S2.4.
√
0.2∆Amin bM/mc
in (42) is non-empty.
For any M ≥ m it holds that Ωm,M
min
Proof. By Lemma S2.3 there exists ω 0 ∈ Ω0.2∆A
. In particular, in holds that
m,m
0
0
min (ASB(ω ), W SB(ω )) ≥ 0.2∆Amin . Define

ω = (ω 0 , . . . , ω 0 , em , . . . , em ) ∈ Ωm,M .
| {z }
bM/mc×

Then ASB(ω), W SB(ω) ≥ 0.2∆Amin

p
bM/mc.

Lemma S2.5.
If Ωδm,M is non-empty for some δ > 0, then there exists ω ∈ Ωδm,M with ASB(ω) = δ.

6

Proof. Fix some ω ∈ Ωδm,M and for 0 ≤  ≤ 1 define ω  ∈ Ωm,M as


ωij

ωij =  ω1j


ωmj + (1 − )ω1j

if i 6∈ {1, m},
if i = 1,
if i = m.

Then W SB(ω  ) ≥ W SB(ω) ≥ δ for all 0 ≤  ≤ 1 and  7→ ASB(ω  ) is continuous
?
with ASB(ω 0 ) = 0 and ASB(ω 1 ) ≥ δ. Thus, there exists an ? ∈ (0, 1] such that ω 
has the desired properties.
Lemma S2.6.
If Ωδm,M in (42) is non-empty, then there exists a quadratic mixing matrix ω ∈ Ωm,m
√
√
such that ASB(ω) = δ∆Amin /(9 M ak ) and W SB(ω) ≥ δ∆Amin /(9 M ak ).
√

Proof. If Ωδm,M is non-empty, it follows from Lemma S2.1 that δ ≤ √M (1+mak ) and,
2m(m+1)
√
hence, δ∆Amin /(9 M ak ) ≤ 0.2∆Amin . Thus, by Lemma S2.3 it follows that
min /(9
Ωδ∆A
m,m

√

M ak )

min
⊃ Ω0.2∆A
6= ∅.
m,m

√
I.e., there exists ω ∈ Ωm,m such that ASB(ω), W SB(ω) ≥ δ∆Amin /(9 M ak ). Now
the assertion follows from Lemma S2.5.

S3 Proof of the lower bound in Theorem 3.1
Proof of Theorem 3.1. The proof of Theorem 3.1 is divided into two steps, corresponding to the two different estimation errors of Π̂ and ω̂, respectively. We start with the
first term on the r.h.s. of the assertion which corresponds to the estimation error of
ω̂. The idea is to construct a hyperrectangle of maximal size which is a subset of
N δ,Λ and then apply results of Donoho et al. (1990) (for fixed selection matrix Π). In
the following, ω ? will denote the center of this hyperrectangle and the matrix E will
?
?
denote the perturbation (of maximal
p size) around ω . ? To this end, let
p ω ∈ Ωm,M
?
be such that ASB(ω ) = 0.2∆Amin bM/mc, W SB(ω ) ≥ 0.2∆Amin bM/mc, and
?
ωmj
≥ 0.4∆Amin /(1 + mak ) for j = 1, . . . , M (existence follows from Lemma S2.4 and
S2.5). For  ∈ (0, 1)(m−1)×M define


11
...
1M


..


.
ω  := ω ? + 
(47)
 = ω ? + E.
 (m−1)1
...
(m−1)M 
Pm−1
Pm−1
− i=1 i1 . . . − i=1 iM
Let  := maxij |ij |. If
0.4∆Amin
≥ (m − 1),
1 + mak

7

(48)

then all entries of ω  are non-negative and ω  ∈ Ωm,M . For ω  to be an element of
Ωδm,M we further need that W SB(ω  ), ASB(ω  ) ≥ δ. To this end, note that

?
kωi· k − kωi−1·
k = kωi·? + Ei· k − kωi−1·
+ Ei−1· k
?
≥ kωi·? k − kωi−1·
k − kEi· k − kEi−1· k
p
√
2
≥
0.2∆Amin bM/mc − m M 
1 + mak

and thus
W SB(ω  ) ≥ 0.2∆Amin

p

bM/mc −

1 + mak √
m M .
2

(49)

Further, note that for e ∈ ∆Am and E as in (47)
2

keEk =

M
X

m−1
X

j=1

i=1

!2
(ei − em )ij

2

≤ M ((m − 1)2ak )

and thus,
ASB(ω  ) = minm ke(ω ? + E)k
e∈∆A

≥ min keω ? k − keEk = ASB(ω ? ) − keEk
e∈∆A
p
√
≥ 0.2∆Amin bM/mc − 2(m − 1)ak  M .

(50)

Summing up, (49) and (50) yield that ω  ∈ Ωδm,M if
p
√
δ ≤ 0.2∆Amin bM/mc − m2 ak M 
(51)
p
√
and (48) holds. As δ ≤ ∆Amin bM/mc(45 2)−1 and M ≥ 2m, (51) and (48) hold
(m−1)×M
for all  ∈ R+
with
∆Amin
≤ √
=: ? .
(52)
6 2m5/2 ak
Now let Π? be the selection matrix such that
Π? A = e 1

e2

em

...

e1

e2

...

>

,

where er ∈ Rm is the r-th unit vector in Rm . Then, as Λ ≤ bn/mc,

Θ := Π? Aω  :  ∈ [0, ? ]m−1×M ⊂ N δ,Λ
and for Π? Aω  ∈ Θ one observes in (14)
i.i.d.

Y1 , . . . , Ybn/mc ∼ N (ω  , σ 2 IM m×M m ).
Define

Θ̃ := ω  :  ∈ [0, ? ]m−1×M .

8

Θ̃ is almost an hyperractangle. To make it a proper hyperractangle, we have to remove
the last column of the matrices in Θ, namely

Θ̃0 := {(ωij
)1≤i≤m−1 :  ∈ [0, ? ]m−1×M }.
1≤j≤M

Note that





inf sup Eθ kθ̂ − θk22 ≥ inf sup Eθ kθ̂ − θk22 .
θ̂ θ∈Θ̃

θ̂ θ∈Θ̃0

Then it follows from (Donoho et al., 1990, (2.1), (3.4) and Proposition 3) that


inf
sup
EΠAω kθ̂ − ΠAωk2
θ̂ ΠAω∈N δ,Λ



sup EΠAω kθ̂ − ΠAωk2
θ̂ ΠAω∈Θ


= bn/mc inf sup Eω kθ̂ − ωk2
≥ inf

θ̂ ω∈Θ̃



≥ bn/mc inf sup Eω kθ̂ − ωk2
θ̂ ω∈Θ̃0

(? )2 σ 2 /bn/mc
(? )2 + σ 2 /bn/mc

−1
1
4
≥ 0.4M
+
,
σ 2 (m − 1) n(? )2

≥ (1.25)−1 bn/mcM (m − 1)

where for the last inequality we used n ≥ 2m. Together with (52) this gives
inf

sup



EΠAω kθ̂ − ΠAωk

2




≥0.4M

θ̂ ΠAω∈N δ

288m5 a2k
1
+
(m − 1)σ 2
(∆Amin )2 n

−1
.

Now we show the second part of the proof which corresponds to the estimation error of Π. The idea is to fix a suitable mixing matrix ω ∈ Ωδm,M and thus, reduce the estimation problem to a classification problem on the finite set of possible selection matrices Π. This can be considered as a testing problem which
√ al0.2
2
lows to apply the Neyman-Pearson lemma. As δ ≤ 9a
(∆A
)
bM/mc/
M it
min
p
p k
follows that 9δak /∆A
M/bM/mc
≤
0.2∆A
bM/mc.
Thus,
by
Lemma
min
√ min
9δak /∆Amin

M/bM/mc

S2.4, Ωm,M
is non-empty and hence, by Lemma S2.6 there exists
a quadratic mixing matrix ω δ ∈ Ωm,m such that
p
p
√
ASB(ω δ ) = (9δak /∆Amin M/bM/mc)/(9 M ak /∆Amin ) = δ/ bM/mc
p
and W SB(ω δ ) ≥ δ/ bM/mc. Hence,
Θ := {ΠA(ω δ , . . . , ω δ , em , . . . , em ) : Π Λ-separable} ⊂ N δ,Λ .
| {z }
bM/mc×

9

(53)

Then for ΠAω ∈ Θ one observes in (14)
i.i.d.

Y1 , . . . , YbM/mc ∼ N (ΠAω δ , σ 2 Inm×nm ).

(54)

For the finite parameter space Θ̃ := {ΠAω δ : Π Λ-separable} it holds that
min kθ − θ0 k2 = min
θ6=θ 0 ∈Θ̃

θ6=θ 0 ∈Θ̃

=

n
X

0 2
kθj· − θj·
k

j=1

min

θ6=θ 0 ∈Θ̃
∃!j ? :θj ? 6=θj0 ?

kθj ? · − θj0 ? · k2 = ASB(ω δ )2 =

δ2
.
bM/mc

Lemma S6.1 yields for any estimator θ̂




sup EΠAω kθ̂ − ΠAωk2 ≥ sup Eθ kθ̂ − θk2
ΠAω∈N δ,Λ

θ∈Θ



= supbM/mc Eθ kθ̂ − θk2

(55)

θ∈Θ̃



≥ δ 2 sup Pθ θ̂ 6= θ .
θ∈Θ̃
2

δ
Now let θ, θ0 ∈ Θ̃ be fixed such that kθ − θ0 k2 = bM/mc
. Then the Neyman-Pearson
PbM/mc
lemma yields for Ȳ := i=1
mYi /M with Yi as in (54) and normally distributed
Z ∼ N (0, σ 2 /bM/mcInm×nm ) that




 1 

Pθ θ̂ 6= θ + Pθ0 θ̂ 6= θ0
sup Pθ θ̂ 6= θ ≥
2
θ∈Θ̃



1 
Pθ θ̂ 6= θ + Pθ0 θ̂ = θ
≥
2

1
≥ inf Pθ kȲ − θk2 − kȲ − θ0 k2 > u
2 u∈R

+ Pθ0 kȲ − θk2 − kȲ − θ0 k2 < u

1
= inf P 2Z > (θ0 − θ) > u + kθ − θ0 k2
2 u∈R

+ P 2Z > (θ0 − θ) < u − kθ − θ0 k2


kθ − θ0 k2
> 0
= P Z (θ − θ) >
2
!
p
 
kθ − θ0 k bM/mc
δ
σ − δ22
=1−Ψ
=1−Ψ
≥
e 8σ ,
2σ
2σ
2δ
where Ψ denotes the cumulative distribution function of the standard normal
and the
√
last inequality follows from Mill’s ratio and 1 − 4σ 2 /δ 2 ≥ 1/2 as δ ≥ σ 8. With (55)

10

this gives



δ2
σ
EΠAω kθ̂ − ΠAωk2 ≥ δ e− 8σ2
2
ΠAω∈N δ,Λ
sup

This finishes the proof.

S4 Proof of the upper bound in Theorem 3.3
Proof of Theorem 3.3. By Theorem 2.2 we can write θ̂ = Π̂Aω̂ in a unique way. We
have that for any ΠAω ∈ N δ,Λ


EΠAω kθ̂ − ΠAωk2




(56)
= EΠAω kθ̂ − ΠAωk2 1{Π̂=Π} + EΠAω kθ̂ − ΠAωk2 1{Π̂6=Π} .
We start with the second
 term. The idea is to bound it with the classification error,
that is PΠAω Π̂ 6= Π , and then apply exact recovery as in Theorem 2.2. As the
entries of the n × M matrices θ̂ and ΠAω are contained in the range of the alphabet
[0, ak ], it follows that




EΠAω kθ̂ − ΠAωk2 1{Π̂6=Π} ≤ a2k nM PΠAω Π̂ 6= Π


X
= a2k nM
PΠAω Π̂ = Π .
Π6=Π
δ
For a fixed
−
√ Λ-separable Π and any ω ∈ Ω it follows from Theorem 2.2 that kΠAω
ΠAωk ≥ Λδ/(1 + mak ) =: c. Further, as by separability rank(ΠA) = m and Ωδ ⊂
Rm×M , there exists a rotation matrix R such that for Θ := {ΠAω − ΠAω : ω ∈
Ωδm,M } ⊂ RnM and Θ̃ := RΘ it holds for all θ ∈ Θ̃ that θmM +2 = . . . = θnM = 0.

11

This gives
PΠAω




Π̂ = Π ≤ PΠAω

!
kY − ΠAωk2 >

min kY − ΠAωk2

ω∈Ωδm,M



= PΠAω kZk2 > min kZ + θk2
θ∈Θ


= P kZk2 > min kZ + θk2
θ∈Θ̃


θ
kθk
= P max −2Z >
−
>0
kθkσ
σ
θ∈Θ̃
!
PmM +1
|θ|
Zi
c
i=1
max
>
≤P
max
i=1,...,M m+1 σ
kθk
2σ
θ∈Θ̃


Z1 √
c
≤ (M m + 1)P
mM + 1 >
σ
2σ


c
= (M m + 1)P |N (0, 1)| > √
2 mM + 1σ
√
2 mM + 1σ − 8(mMc2+1)σ2
≤ (M m + 1)
,
e
c
where we considered the noise matrix Z ∈ Rn×M in (14) as a vector Z ∈ RnM (with
i.i.d.

entries Zi ∼ N (0, σ 2 ), i = 1, . . . , nM ) and for the last inequality we used Mill’s
ratio. As the number of Λ-separable selection matrices Π is bounded by nmk m , it
follows that


Λδ 2
(1 + mM )3/2 − 8(mM +1)(1+ma
2 2
k) σ .
√
PΠAω Π̂ 6= Π ≤ 2σnk m m(1 + mak )
e
Λδ
and


Λδ 2
M 5/2 −
EΠAω kθ̂ − ΠAωk2 1{Π̂6=Π} ≤ 12σn2 k m m7/2 a3k √ e 8(mM +1)(1+mak )2 σ2 .
Λδ
This gives the second term of the r.h.s. of the assertion.
Now we consider the first summand on the r.h.s. of (56). The idea is to bound
the minimax risk conditioned on Π̂ = Π with the minimax risk of the LSE on the
linear subvector space imag(ΠA). To this end, let N δ (Π) ⊂ N δ,Λ denote the set of
all ΠAω ∈ N δ,Λ with Π = Π. Further, let θ̂0 ∈ argminΠAω∈N δ (Π) kY − ΠAωk2 be the
least-squares estimator restricted to Π̂ = Π. Then, clearly, θ̂ = θ̂0 on {Π̂ = Π} and
thus




EΠAω kθ̂ − ΠAωk2 1{Π̂=Π} = EΠAω kθ̂0 − ΠAωk2 1{Π̂=Π}


≤ EΠAω kθ̂0 − ΠAωk2 .

12

Thus, for all fixed Π as in (15)
sup
ΠAω∈N δ (Π)

≤

sup



EΠAω kθ̂ − ΠAωk2 1{Π̂=Π}


EΠAω kθ̂0 − ΠAωk2 .

ΠAω∈N δ (Π)

Clearly, N δ (Π) ⊂ imag(ΠA)M with dim(imag(ΠA)M ) = mM . Thus, for the LS
estimator on imag(ΠA)M , θ̂00 ∈ argminθ̃∈imag(ΠA)M kY − θ̃k2 , it follows from Lemma
S6.2 that


sup
EΠAω kθ̂0 − ΠAωk2
ΠAω∈N δ (Π)

≤4

sup



Eθ kθ̂00 − θk2 = 4σ 2 mM,

θ∈imag(ΠA)M

which finishes the proof.

S5 Proofs for the results on the estimation error in
Section 4
Proof of Theorem 4.2. The first√assertion follows directly from the first part of Lemma
S6.3 with  = kΠAω−Π0 Aω 0 k/( nmak ). The second assertion follows from the second
part of Lemma S6.3 with  % kΠAω − Π0 Aω 0 k.
Proof of Theorem 4.3. It follows directly from combining the first part of Theorem 4.2
and Corollary 3.4 that
 
2 
c2
√
c1
1
− 81 σM2
inf sup EΠAω d (Π, ω), (Π̂, ω̂)
M
.
& σ2 M
+
σ
e
nma2k
nm2 a2k
Π̂,ω̂ Π,ω
n
o
√
For B := kΠAω − Π̂Aω̂k ≤ c1 M /(1 + mak ) , Markov’s inequality gives


√
PΠAω (B c ) = PΠAω kΠAω − Π̂Aω̂k > c1 M /(1 + mak )


EΠAω kΠAω − Π̂Aω̂k2
.
≤
c21 M/(1 + mak )2
√
√
As for any ΠAω, ΠAω ∈ N δ,Λ , d ((Π, ω), (Π, ω)) ≤ M (1/m+1) ≤ M 3/2 it follows

13

from combining the second part of Theorem 4.2 and Corollary 3.4 that
 
2 
inf sup EΠAω d (Π, ω), (Π̂, ω̂)
Π̂,ω̂ Π,ω

 

2
≤ inf sup EΠAω d (Π, ω), (Π̂, ω̂) 1B + 9/4M PΠAω (B c )
Π̂,ω̂ Π,ω





≤ inf sup EΠAω kΠAω − Π̂Aω̂k2 + 9/4M



EΠAω kΠAω − Π̂Aω̂k2

c21 M/(1 + mak )2


9/4(1 + mak )2
2
= inf sup EΠAω kΠAω − Π̂Aω̂k
1+
c21
Π̂,ω̂ ω,ω


 ma 2

k
. inf sup EΠAω kΠAω − Π̂Aω̂k2
c
Π̂,ω̂ ω,ω
1
c2 c2
2 2
3 2
M
√
1
m a −
m ak
+ M σ √ k3 e 16m3 (1+ak )2 σ2 .
. σ2 M
2
c1
c2 c1
Π̂,ω̂ Π,ω



S6 Additional lemmas
Lemma S6.1.
For a finite parameter space Θ ⊂ Rn and any estimator θ̂


supθ∈Θ Eθ kθ̂ − θk2

 ≤ max kθ0 − θ00 k2 .
min kθ0 − θ00 k2 ≤
θ 0 6=θ 00
θ 0 6=θ 00
supθ∈Θ Pθ θ̂ 6= θ)
Proof. It holds that




X
sup Eθ kθ̂ − θk2 = sup
kθ̃ − θk2 Pθ θ̂ = θ̃
θ∈Θ

θ∈Θ

θ̃∈Θ\θ

= sup Pθ (θ̂ 6= θ)
θ∈Θ

X

kθ̃ − θk2

θ̃∈Θ\θ



Pθ θ̂ = θ̃
Pθ (θ̂ 6= θ)

,

where, for every θ ∈ Θ
min kθ0 − θ00 k2 ≤

X

θ 0 6=θ 00

kθ̃ − θk2



Pθ θ̂ = θ̃

θ̃∈Θ\θ

14

Pθ (θ̂ 6= θ)

≤ max
kθ0 − θ00 k2 .
0
00
θ 6=θ

Lemma S6.2.
Let V be a subvector space of Rd , A ⊂ V an arbitrary subset, and Y ∈ Rd . Further,
let θ̂A (Y ) ∈ argminθ̃∈A kY − θ̃k and θ̂V (Y ) ∈ argminθ̃∈V kY − θ̃k. Then
kθ̂A − θk2 ≤ 4 kθ̂V − θk2 .

∀θ ∈ A :

(57)

Proof. Let θ ∈ A be fixed. If kθ̂A − θk2 = 0, (57) holds trivially. Further, if Y ∈ A,
it holds that θ̂V = θ̂A = Y , and hence, (57) follows trivially, too. So assume that
kθ̂A − θk2 > 0 and Y 6∈ A.
Choosing an appropriate coordinate system, we may w.l.o.g. assume that
V = {x ∈ Rd : x1 = . . . = xr = 0},
with dim(V ) = d − r. Let pr be the orthogonal projection onto V , i.e.,
pr : (x1 , . . . , xd )> 7→ (xr+1 , . . . , xd )> .
Then,
argminθ̃∈A kY − θ̃k2 = argminθ̃∈A

r
X

Yi2 + k pr(Y ) − pr(θ̃)k2

i=1

= argminθ̃∈A k pr(Y ) − pr(θ̃)k2
and analog argminθ̃∈V kY − θ̃k2 = argminθ̃∈V k pr(Y )−pr(θ̃)k2 , kθ̂V −θk2 = k pr(θ̂V )−
pr(θ)k2 , and kθ̂V − θk2 = k pr(θ̂V ) − pr(θ)k2 . Thus, we may w.l.o.g. assume that
V = Rd , i.e., θ̂V = Y . Then
kY − θk2
kθ̂A − θk2

≥

min

kx − θk

!2
≥

kθ̂A (x) − θk
!−2
1
kx − θ̂A (x)k
≥ ,
= 1 + max
x6∈A
kx − θk
4
x6∈A

min
x6∈A

kx − θk
kθ̂A (x) − xk + kx − θk

where the last inequality follows from the definition of θ̂A .
Lemma S6.3.
Let ΠAω, Π0 Aω 0 ∈ N δ , then for all  > 0
1.
kΠAω − Π0 Aω 0 k ≥

√

nmak 

⇒

d ((Π, ω), (Π0 , ω 0 )) ≥ ,

2. if  < δ/(1 + mak ), then
kΠAω − Π0 Aω 0 k < 

⇒

15

d ((Π, ω), (Π0 , ω 0 )) < .

!2

Proof. From kΠAω − Π0 Aω 0 k ≥

√

nmak  it follows that

max k(ΠAω)j· − (Π0 Aω 0 )j· k ≥ mak 

j=1,...,n

and
 ≤ max k(ΠAω)j· − (Π0 Aω 0 )j· k/(mak ) ≤
j=1,...,n

√
M /m.

Hence, by Theorem 2.2 maxi=1,...,m kωi· − ωi·0 k ≥  or Π 6= Π0 . Thus, it follows
that d ((Π, ω), (Π0 , ω 0 )) ≥ , which shows the first assertion. If  < δ/(1 + mak )
and maxj=1,...,n k(ΠAω)j· − (Π0 Aω 0 )j· k ≤ kΠAω − Π0 Aω 0 k <  it follows from Theorem 2.2 that maxi=1,...,m kωi· − ωi·0 k <  and Π = Π0 and thus it follows that
d ((Π, ω), (Π0 , ω 0 )) < , which shows the second assertion.

16

