I have a system whose output I am analyzing. I've used a very simple hypothesis test with great results and now am curious why did I get such great results -- did I perhaps stumble on some kind of an "optimal" test without much formal statistics training?

My system can be in two states. When it's in state 0 the output can is independently drawn from standard normal distribution $\mathcal{N}(0,1)$ (well, technically, very, very close approximation to it-- but it's normal enough for my purposes). When it's in state 1, the output is still independently normally-distributed with variance 1 but the means are perturbed. Thus, the output $i$ in state 1 has a mean-shifted normal distribution $\mathcal{N}(\mu_i,1)$. For large number $n$ of observations of the system in state 1, the average $\frac{1}{n}\sum_{i=1}^n\mu_i=0$. I also know that quantity $\frac{1}{n}\sum_{i=1}^n\mu_i^2=M<\infty$ since this is a real-world system (though I do not know the exact value of $M$, just the fact that it's finite). Besides those two facts, there is no structure to means in state 1.

I devised the following hypothesis test to determine whether the system is in state 0 or 1, where the null hypothesis $H_0$ corresponds to the system being in state 0, and alternate $H_1$ to state 1. First, I collect a sequence of $n$ observations $\{x_i\}_{i=1}^n$. Then I compute the following test statistic:

$$S_n=\left(\frac{1}{n}\sum x_i^2\right)-1$$

I then pick a threshold $t>0$ and accept null hypothesis if $S_n<t$, rejecting it if $S_n\geq t$.

This test just made sense to me and is surprisingly easy to analyze (at least surprising to me without much formal background in mathematical statistics). The mean of $S_n$ when the system is in state 0 is 0, and the variance is $2/n$. In state 1 the mean of $S_n$ is $M$ and the variance is $(4M+2)/n$. Thus, I can upper-bound both the false-positive error probability $\alpha$ and the probability that I will accept the null hypothesis in error $\beta$ for a given $n$ and $M$ using the Chebyshev's Inequality:

$$\begin{array}{rcl}\alpha&\leq&\frac{2}{nt^2}\\ \beta&\leq&\frac{4M+2}{n(M-t)^2} \end{array}$$

Using the two inequalities I can pick an appropriate threshold $t$ and the number of observations $n$ I need to classify the state of the system for some particular value of $M$ given some tolerance for errors $\alpha$ and $\beta$. The math is simple and completely makes sense: for some fixed $\alpha$, as I increase $n$, my threshold $t$ decreases, and so does $\beta$.

Is the reason for this beauty and simplicity that the test I devised is in some way "optimal"? If so, can the community suggest a way prove this optimality? I did find "Uniformly Most Powerful tests" in a stats book, but I am not sure if I have that (the book was quite confusing)...

  • 2
    $\begingroup$ I think you have discovered a form of what is known in the engineering literature as a square-law detector or noncoherent detector. If the sequence of $\mu_i$ is known, the optimum detector uses the statistic $\sum_ix_i\mu_i$, heuristically $\mu_i$ and $x_i$ are more likely to have the same sign. When the $\mu_i$ are not known, the detector essentially makes an "energy measurement" $\sum_i x_i^2$ (you might remember $V^2/R$ from high-school or college physics) with large energy meaning signal is present and little energy meaning signal is absent. $\endgroup$ Dec 21 '11 at 2:30
  • $\begingroup$ Dilip, interesting point about energy measurement, yes, I remember the $P=V^2/R$ formula from high school physics. I've posted on stats.SE before, though I thought this question is general enough to be posted here. How can I ask the moderators to move this question to stats.SE? $\endgroup$
    – M.B.M.
    Dec 21 '11 at 18:48
  • $\begingroup$ How do you compute a variance of $2/n$ in state 0 and $(4M+2)/n$ in state 1? (I obtain different values for the variances.) Also, how can you be sure $\frac{1}{n}\sum_{i=1}^n\mu_i$ is exactly $0$ and not, say, $O(1/n)$ or $o(1/n)$? How do you know the values of the $\mu_i$ (and, thence, the value of $M$)? Unless you know them exactly, independently of the data, Chebyshev's inequality does not apply here. And if you do know the $\mu_i$ independently of the data, why do you need any kind of test at all? $\endgroup$
    – whuber
    Dec 22 '11 at 3:28
  • $\begingroup$ Re: $\mu_i$'s -- they are not random, but they are unknown. My uncertainty about $S_n$ comes only from the underlying independent normal distribution. The variance of $S_n$ in state 0 is its 4th moment minus the second moment squared. For each squared observation it's 3-1=2. To compute the variance in state 1, I first find the variance of an individual squared observation (4th moment of mean-shifted normal with variance 1) then subtract the squared-mean of the squared observation. Since the individual squared observations are independent, variance is additive... $\endgroup$
    – M.B.M.
    Dec 22 '11 at 5:03
  • $\begingroup$ @DilipSarwate: can you provide a reference on why $\sum_i x_i \mu_i$ is the optimal test? I've encountered a similar problem in another context. Thanks! $\endgroup$
    – J Li
    Feb 13 '16 at 17:20

You can do simpler than using Chebyshev’s inequality : the distribution of $\sum_{i=1}^n X_i^2$ under $H_0$ known, it is a $\chi^2(n)$. From this you can get a threshold $t_\alpha$, using simply a chi square table: $t_\alpha = x_{n,1-\alpha}$, the quantile of level $1-\alpha$ of $\chi^2(n)$.

This is a particularly simple situation as the variance of the law under $H_0$ in known. This is almost like the kind of artificial stuff I teach to beginners, I would be pleased if you were kind enough to describe the concrete situation that lead you to this model.

You can compare the threshold you get from Chebyshev inequality: for $\sum X_i^2$ this $n\left(1+\sqrt{{2n\over\alpha}}\right)$ to the $t_\alpha$ obtained by quantiles of $\chi^2(n)$, you will see that the quantiles are much lower; this will be much more powerful than you thought.

Post Scriptum To adress your original question: It is not easy to show that this test is optimal in any sense; the usual tool is Neyman-Pearson lemma but it is for tests where both $H_0$ and $H_1$ are point hypotheses, which is not the case there. Other tools are avalaible when the parameter space is $\mathbb{R}$ but here the parameters are $\mu_1, \dots, \mu_n$... So this is beyond my abilities. A family of test statistics that could be investigated is $T_\lambda = \sum_i |X_i|^\lambda$. Intuitively, I would say that the value of $\lambda$ giving the most powerful test depends on the values $(\mu_i)_{i=1,\dots,n}$. If you knew these values, the uniformely most powerful test would be given by the likelihood ratio test, this is Neyman-Pearson Lemma.


Correct me if I'm wrong because I'm still learning this material. The test above should be most powerful for some level $\alpha$.

The distribution under the alternative should be non-central Chi-Squared with $n$ degrees of freedom and non-centrality parameter $nM$, $\chi^2_n(nM)$ because if the system is in state $i$ then the entire sequence is generated in that state and all we want to know is which state the system is in.

Given this, and the fact that the distribution under the null is $\chi^2_n$, both hypotheses are simple and the Neyman-Pearson Lemma can be applied. According to Neyman-Pearson, a most powerful test is given by: reject if $p_1/p_0 \geq c$ and do not reject if $p_1/p_0 < c$ (where $p_0$ and $p_1$ are the pdfs of $\sum{X_i^2}$ under the null and alternative respectively).

The ratio of the pdfs is $r(y) = e^{-nM/2} {_0F_1}(;\frac{n}{2};\frac{nMy}{4})$ (where ${_0F_1}$ is a generalized hypergeometric function, and $y = \sum{x_i^2}$).

To find a most powerful test of size $\alpha$, we could choose $c$ such that $\mathbb{E}_{\chi^{2}_n}[I_{r(y)<c}] = \alpha$. But because $r(y)$ increases with $y$, this is equivalent to finding $c'$ such that $\mathbb{P}_{\chi^2_n}(y<c') = \alpha$.

$S_n$ is a normalized and shifted version of $\sum{x_i^2}$, applying the same transformations to $c'$ gives $t$ such that rejecting when $S_n \geq t$ and not rejecting when $S_n < t$ gives a most powerful test at level $\alpha$. The bound derived above works for bounding the power, but since we know the distribution under the alternative we can figure out the power of the test exactly.

  • $\begingroup$ Thanks for a nice answer. I think that what I was wondering about in the original question is why the square-law detector (the term Dilip Sarwate used to describe my test in his comment) works so well. Besides the test statistic having well-studied distributions under either hypothesis, can we say something about the test itself? While obviously, using the mean of the observation yields no information, why can't we do better by taking, each observation to, say the power of 1.2 (using the magnitudes of complex numbers) or, say power of 4? Why does squaring work so well? $\endgroup$
    – M.B.M.
    Jan 4 '14 at 19:14
  • $\begingroup$ Very welcome. No other power would work as well because the decisions by the MP test is unique by the NP lemma. What intuition I can find is, consider the distribution of the data in n-dimensions (say 3). The null distribution is normal centered at 0 and the alternative could be centered anywhere on some circle of Euclidean radius nM. The regions that the square norm is constant on are spheres, which is only true for the square norm. Since the shape of these regions determines reject-no reject decisions, we want the most symmetric shape because the alternative is centered on a circle. $\endgroup$ Jan 5 '14 at 8:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.