Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
118,997 characters · 29 sections · 31 citation commands
Improved Tests for Mediation
\noindentKeywords: similar test, augmented LR test, coherence.
\spacingset{1.9}
Testing for a mediation effect has important applications in many disciplines, including psychology, sociology, epidemiology, accounting, marketing, economics and business.\footnote{ See, for example, Baron1986, Coletti2005, Mackenzie1986, Alwin1975, Freedman1992, Heckman2015a, Heckman2015b.} A simple context for the problem - which we will use to motivate the results to follow - is a model of the type
This model is generally given a causal interpretation with $x$ exerting a causal influence on $y$ via the mediating variable $m$. The influence of $x$ on $y$ may be both direct (the term $\tau x),$ and/or indirect via the term $ \theta _{2}m$ if $\theta _{1}$ is nonzero in the first equation. There is no mediation effect if either $\theta _{1}=0$ such that $x$ does not appear in ( (ref)), or $\theta _{2}=0$ such that $m$ does not appear in ((ref)). So a test for the absence of a mediation effect is a test of the composite null hypothesis $H_{0}:\theta _{1}\theta _{2}=0$. Controls can be added to the model to avoid unmeasured confounding effects. The vectors $y,x,$ and $m$ are then the residuals after regression on these controls. This has no bearing on what follows if standard assumptions regarding controls are met: after including controls (covariates) no un-measured common causes exist for the relations between (i) $x$ and $m,$ (ii) $x$ and $y,$ (iii) $m$ and $y$, and (iv) $m$ and $y$ that are affected by $x$; see for instance vdweele2015.
This testing problem is complicated by the fact that, even asymptotically, there is a nuisance parameter present under the null - either $\theta _{1}$ or $\theta _{2}$ may be nonzero - and this seriously impacts the properties of most of the tests that have been proposed for the problem; see Mackinnon2002 for a survey. Typically, the extant tests have very poor power behavior near the origin ($\theta _{1}=\theta _{2}=0).$ Specifically, the null rejection probability (NRP) of the test can be very much smaller than its nominal size - near zero in fact - and its power and NRP can be very nearly equal. The bank of standard tests exhibiting this behavior all reject the null hypothesis when some particular test statistic is large. However, in a recent paper, VG2-2021, VG2 hereafter, have shown that both NRP and power can be improved considerably - particularly near the origin - by using a critical region that cannot be defined in this way, but is simply a subset of a two-dimensional sample space. After reducing the problem by invariance, they consider a critical region $(CR)$ consisting of the likelihood ratio region ($CR_{LR})$, augmented by an additional region closer to the origin. This additional region is carefully constructed using a piecewise-linear spline, and is optimized in terms of both NRP and power.
In this paper we employ the same idea - augmenting the $CR_{LR}$ by an additional region - but, in the interests of pragmatism, our focus here will be on constructing a test that is very simple to apply in empirical research, yet has\ NRP very nearly constant.
To motivate the proposed test we first show that, for certain test sizes (including the popular choices $\alpha =.01,.05,$ and $.10)$ a test exists that is exactly similar asymptotically, i.e. with constant NRP equal to $ \alpha $, and we show how to construct it. The construction of this exact test resembles those mentioned by lehmann1952testing and later by nomakuchi1987note; see also berger1989 for related constructions. As in those earlier examples, however, the critical regions have some undesirable characteristics including the lack of coherency. A test procedure is said to be coherent if, when the test rejects the null at level $\alpha $, it also rejects at all levels greater than $\alpha $. The likelihood ratio test has this property, see VG2-2022, but we show below that the exact test does not.
We therefore propose an alternative, easily constructed, test that is close to being exact, and which avoids some of the undesirable aspects of the exact test. We call this the "simply-augmented LR test" as it adds a region with linear boundary defined by two points to the $CR_{LR}$. Specifically, using the two common $t$-statistics for $\theta _{1}$ and $\theta _{2},$ if $ v_{1}=\min \{t_{1}^{2},t_{2}^{2}\}$ and $v_{2}=\max \{t_{1}^{2},t_{2}^{2}\},$ then the $.05$ level test is simply:\newline reject $H_{0}:\theta _{1}\theta _{2}=0$ if $v_{1}>$ $ 3.841$ (the LR test) or $v_{1}/v_{2}>0.8744$ (augmentation), \ \newline i.e. reject if a traditional LR\ test rejects and otherwise check the augmented region. Only one number is required - in addition to the usual $ \chi _{1}^{2}$-critical value - in order to implement the test. We provide these critical values for every percentile level in Table (ref) at the end of the paper. We prove theoretically that this test cannot be exactly size-correct. It can be modified to ensure that its NRP $\leq \alpha $ for all parameter values under the null, but the resulting (truncated-augmented) test lacks coherency. Coherency is crucial for the interpretability and derivation of p-values in empirical research. Table (ref) also shows the straightforward calculation of p-values based on the coherent simply-augmented LR test. The new test is far superior to the LR test in terms of NRP, and, trivially (because their critical regions are larger), also has greater power.
For any test with critical region $w$, and where the distribution of the statistics involved depends on a vector of parameters $\psi ,$ we denote the power of the test by $P_{w}(\psi ),$ and the NRP when the null distribution depends on the parameter $\psi _{0}$ by $P_{w}(\psi _{0}).$ The size of the test is as usual defined to be $\sup_{\psi _{0}}P_{w}(\psi _{0}).$ We emphasize that the issue we are concerned with here is not that of finding tests of the correct size - the LR and Sobel's Wald test both have this property - but that the usual tests can have NRP and power that are near zero in relevant parts of the parameter space where mediation effect is small or imprecisely estimated.
Section (ref) motivates restricting our search for an improved testing procedure to the squared $t$-statistics. We formulate the asymptotic problem under minimal assumptions that allow for\ nonnormality and heteroskedasticity, which is important empirically. Section (ref) analyzes the LR test and shows the poor NRP and power properties, motivating the search for an improved test. Section (ref) proves the existence of an exact, but incoherent test. Section (ref) considers the augmented tests and derives their properties. Section (ref) shows the coherency of the simply-augmented test and how this leads to p-values. Section (ref) shows the relevance of the (asymptotic-based) procedures in finite samples encountered in practice by simulation and an empirical application before concluding in Section (ref).
The testing problem for $H_{0}:\theta _{1}\theta _{2}=0$ respects several symmetries, such as the signs of the coefficients. Furthermore, the truth or falsity of the null is irrespective of the error variances, or the value of $ \tau $. We are looking for testing procedures that are invariant to these transformations. They can be based solely on the ordinary $t$-ratios for $ \theta _{1}$ and $\theta _{2}$ since they are shown to be the maximal invariants under the relevant group of transformations in Theorem A.1 of the Appendix under Gaussianity. The LR test and Sobel's Wald test are in fact basic functions of only these two statistics. The Gaussian assumption is unnecessarily restrictive however. The two $t$-statistics converge to normality under the much weaker conditions based on White1980, that we introduce next. Heteroskedasticity in particular is ubiquitous in empirical research, and can be addressed using robust $t$-statistics. Proposition (ref) shows the convergence of the two relevant robust $t$-statistics defined in Equation ((ref)) below to a normal distribution under the following assumptions.
The essential conditions are that the model defined by the regressions ((ref)) and ((ref)) are correct and therefore that $u_{1}$ conditional on $x$ and $u_{2}$ conditional on $\left( x,m\right) $ have expectations zero. By the law of iterated expectations this implies that the disturbances $u_{1}$ and $u_{2}$ are uncorrelated, as is usually assumed explicitly. Observations are assumed to be independent and sufficient higher-order moments should exist. Note that one of the variables in $\mathbf{X}$ is $m$ which satisfies Equation ((ref)). Assumptions can therefore be imposed on $x$ and $u_{1},$ rather than $m$ itself. Further note that Assumption (ref) (iv) implies that the elements of $ \Omega _{n}$ are uniformly bounded and together with (v) ensure uniform boundedness of $\Omega _{n}^{-1}$ (see White1980). Similarly for $Q_{n}^{-1}.$
The remainder of the paper will be based on asymptotic distribution ((ref)) of the robust $t$-statistics. The problem then becomes: we observe independent random variables $T_{1},T_{2},$ with $T_{i}\sim N(\mu _{i},1),i=1,2,$ and wish to test the hypothesis $H_{0}:\mu _{1}\mu _{2}=0.$ It is clear that this problem is invariant under the group of sign changes \thinspace $T_{i}\mapsto -T_{i}$, $i=1,2,$ and under this group of transformations the statistics $f_{i}=T_{i}^{2}$, $i=1,2,$ are maximal invariants. These are independent noncentral $\chi _{1}^{2}$ variates with noncentrality parameters $\lambda _{i}=\mu _{i}^{2}$, $i=1,2.$ Therefore the no-mediation restriction $\theta _{1}\theta _{2}=0$ versus $\theta _{1}\theta _{2}\neq 0$ is equivalent to testing
This problem is clearly also invariant under the group of permutations of $ (f_{1},f_{2}),$ and maximal invariants under this action are $(v_{1},v_{2}),$ with $v_{i}=f_{(i)}$ the $i^{th}$ order statistic (so $v_{2}\geq v_{1}\geq 0).$ Thus, we are finally led to focus attention on the pair of order statistics $(v_{1},v_{2})=(f_{(1)},f_{(2)}),$ which live on the octant $ V=\{(v_{1},v_{2});0\leq v_{1}\leq v_{2}<\infty \}.$ The reader should bear in mind, though, that any test i.e. $CR$ formulated in terms of $ (v_{1},v_{2})$ can equally well be re-expressed in terms of the $t$ -statistics $(T_{1},T_{2})$; see Figure (ref).
The null hypothesis is composite and involves a nuisance parameter $\lambda =\max \{\lambda _{1},\lambda _{2}\},$ the (possibly non-vanishing) noncentrality parameter. It is therefore not obvious how to construct a test (critical region) whose NRP $P_{w}(\lambda )$ does not depend on $\lambda .$ However, we will show below that for each level $\alpha =$ $(r+2)^{-1}$, where $r$ is a non-negative integer, there exists an exact similar test. Trivially, these tests have power functions uniformly above that of the LR test.
As already remarked, the NRP and power of the LR test, and other standard tests, can in fact be extremely small - NRP when the nuisance parameter is small, power when both $\lambda _{1}$, $\lambda _{2}$ are small. The popular Sobel (Wald) test is uniformly much worse than the LR test in both respects and we will therefore not discuss it further.\footnote{ \ The Wald test rejects when $W=\frac{v_{1}v_{2}}{v_{1}+v_{2}}$ is large, using the same critical value as the LR test, but $LR>W.$ The LR is less biased and more powerful. The $W$ power function stays closer to zero longer as $\lambda $ increases. If $\lambda =0$ and $\alpha =.05,$ the Wald NRP is $ .00009$, while $(.05)^{2}=.0025$ for LR test.} There is clearly an incentive to seek a test whose NRP is closer to the nominal size for all values of the nuisance parameter, and has better power. This is the motivation for what follows, which builds on the LR test that we discuss next.
In this section we derive the NRP and power of the likelihood ratio (LR) test. This requires the joint distribution of the maximal invariants, the order statistics $(v_{1},v_{2})$ that play a central role in the derivation and the properties of proposed improved tests. The LR test is derived by minimizing $(t_{1}-\mu _{1})^{2}+(t_{2}-\mu _{2})^{2}$ subject to the constraint $\mu _{1}\mu _{2}=0.$ It is straightforward to show that this results in the following critical region in the space of the order statistics $(v_{1},v_{2})$: reject $H_{0}:\mu _{1}\mu _{2}=0$ when
is large. As usual, the LR test embodies all invariance properties of the testing problem. The critical region for the LR test of nominal size $\alpha $ is given by the set $CR_{LR}=\{\chi _{\alpha }^{2}<v_{1}<v_{2},v_{2}>\chi _{\alpha }^{2}\}$ with $\chi _{\alpha }^{2}$ the $\alpha $ critical value from the $\chi _{1}^{2}$ distribution. Let $g(v;\lambda )$ and $G(v;\lambda ) $ denote the pdf and CDF of the noncentral $\chi _{1}^{2}\left( \lambda \right) $ distribution. The joint distribution of $(v_{1},v_{2})$ as used for all subsequent results in the paper equals:
based on the premise that $T\sim N(\mu ,I_{2})$; see Proposition A.1 in Appendix A and its specialization to the null case. Given this distribution, the following proposition provides a very direct description of the properties of the LR test.
Since, for any fixed $z>0,$ $1-G(z;\lambda )$ is an increasing function of $ \lambda ,$ tending to one as $\lambda \rightarrow \infty $, we arrive at the following corollary:
Thus, the LR test has a size of $\alpha $. However, for small $\lambda ,$ the NRP of the LR test can be as small as $\alpha ^{2}$, and it only approaches the nominal size $\alpha $ as $\lambda \rightarrow \infty $.
The tests we consider below are constructed by augmenting the $LR$ critical region, and it is clear from the expression for $P_{CR_{LR}}(\lambda )$ above that the region added should have null content either exactly equal to $\alpha G(\chi_{\alpha }^{2};\lambda )$ for all $\lambda ,$ rendering the test exact, or have this property approximately. Both exact and approximate augmented LR tests will be constructed below.
The poor NRP and power properties of the classical tests motivate the search for more satisfactory tests. Specifically, we would hope to be able to construct tests whose NRP is $\alpha $, or nearly so, for all $\lambda ,$ and whose power improves on that of the LR test, in particular. In this section we shall show that an exact test does indeed exist for certain choices of $\alpha ,$ and is easily constructed. We confine attention to tests whose critical regions properly contain that of the LR test. That is, if $H_{0}$ is rejected by the LR test it must also be rejected by the new test (but not vice versa).
It is convenient for the derivation of the exact test to introduce a partition of the sample space, the octant $V,$ into three disjoint regions determined by a scalar $z>0$:$\newline $ $\ A_{1}=\{v_{2}>z,z<v_{1}<v_{2}\},\ A_{2}=\{v_{2}>z,0<v_{1}<z\},\ A_{3}=\{v_{2}<z,0<v_{1}<v_{2}\}.\newline $ The first of these, $A_{1},$ is the level-$\alpha $ $CR_{LR}$ when $z=\chi _{\alpha }^{2}$ with acceptance region $AR_{LR}=A_{2}\cup A_{3}$; see Figure (ref) for a graphical comparison between the CR based on the $t$ -statistics $t_{1}$ and $t_{2}$ and the CR derived from the order statistics of the squared $t$-statistics $v_{1}$ and $v_{2}$. In what follows the regions $A_{1}$, $A_{2}$, $A_{3}$ will be assumed to be defined by $z=\chi _{\alpha }^{2}$ with probabilities:
The next proposition determines the NRP of the $CR_{LR},$ i.e. $A_{1}$, augmented by $A_{3}.$
To illustrate the general result, consider first choosing a single value $ z_{1}<\chi _{\alpha }^{2}$, with $\alpha $ to be determined also, and using this to define two disjoint triangular subsets of $A_{3}$: $ A_{30}=\{0<v_{1}<v_{2},0<v_{2}<z_{1}\}$ and $A_{31}= \{z_{1}<v_{1}<v_{2},z_{1}<v_{2}\,<\chi _{\alpha }^{2}\}$. These have combined null probability content
which differs from the target value $\alpha G(\chi _{\alpha }^{2};\lambda )$ by
Since $\alpha =1-G(\chi _{\alpha }^{2}),$ we can choose the pair $ (z_{1},\chi _{\alpha }^{2})$ so that the coefficients of the two noncentral distribution functions both vanish, yielding a test of size $1-G(\chi _{\alpha }^{2})$ for all $\lambda .$ This requirement produces two linear equations, $2G(\chi _{\alpha }^{2})-G(z_{1})=1,$ and $G(\chi _{\alpha }^{2})-2G(z_{1})=0,$ with unique solution $G(z_{1})=1/3,G(\chi _{\alpha }^{2})=2/3,$ so that $\alpha =1-G(\chi _{\alpha }^{2})=1/3.$ This is the case $r=1,$ $\alpha =1/3.$
Generalizing this construction, one may prove, as we do in Appendix A:
When the $z_{i}$ are chosen in this optimal fashion we denote the augmenting region simply by $w_{r}$.
The power function of the exact test described above is obviously higher than that of the LR test for all $(\lambda _{1},\lambda _{2}),$ whatever the value of $r.$ It is easy to check that the probability content of the region $A_{3}$ under the alternative is $G(\chi _{\alpha }^{2};\lambda _{1})G(\chi _{\alpha }^{2};\lambda _{2}).$ The power added by the augmenting region is therefore given by
This is naturally symmetric in $(\lambda _{1},\lambda _{2}),$ and vanishes in the limit as either noncentrality parameter goes to infinity. That is, there is no power gain over the LR test in the limit, but there certainly is for $(\lambda _{1},\lambda _{2})$ near the origin. The NRP gain at the origin is obviously $\alpha (1-\alpha )=(r+1)\alpha ^{2}.$ The power function behaves similarly for points $(\lambda _{1},\lambda _{2})$ close to the origin with a power gain of around $.0475$ over the LR test when $\alpha =.05$.
For a restricted, but relevant, range of nominal sizes ($\alpha $ of the form $(r+2)^{-1})$ the construction described above provides, for the first time, a non-randomized exact test of the no-mediation hypothesis. Whilst the augmenting critical region does contain points close to the origin, which might be considered counter-intuitive, over 90% of the area of the augmenting critical region in the case of $\alpha =.05$ is accounted for by the four largest triangular regions, and these regions are well away from the origin. Nevertheless, the critical region of the exact test does have several undesirable properties. First, the region $CR_{LR}\cup w_{r}=CR_{r},$ say, is not monotone in $(v_{1},v_{2}).$ That is, $(v_{1},v_{2})\in CR_{r}$ does not imply that $(v_{1}^{\prime },v_{2}^{\prime })\in CR_{r}$ when $ v_{1}^{\prime }\geq v_{1}$ and $v_{2}^{\prime }\geq v_{2}.$ Similarly, the acceptance region for the test is not convex, which is somewhat counter-intuitive. Also, unlike the LR test itself, the exact test does not possess an important coherence property, namely, that rejection at level $ \alpha $ does not imply rejection at every level higher than $\alpha .$ That is, as the reader may easily confirm, the critical region for $\alpha =(r+2)^{-1}$ is not a subset of that for $\alpha =(r+1)^{-1}.$ This also rules out the use of p-values. The augmented LR tests introduced in the next section will address some of these deficiencies.
To motivate the class of tests we consider next, observe that one could approximate the exact augmenting region (i.e., the blue $w_{r}$ triangles in Figure (ref)) with the region above a line $ v_{1}=bv_{2} $, for some suitable choice of $b$. For instance, the (geometric) area of the augmenting region of the exact test in the case $ r=18 $ $(\alpha =.05)$ is $1.08403$, which is equal to the area above the line $v_{1}=(.87187)v_{2}$ (red solid line in Figure (ref)) and below $v_{1}=\chi _{.05}^{2}$.\ Unlike the exact CR, the NRP of the region defined in this way does vary slightly with $ \lambda .$ VG2 discusses (in present notation) more general augmenting regions close to the $v_{1}=v_{2}$ line, but bounded below by a piecewise-linear spline. The exact test just described is of this form, but a very special case: the alternate knots are constrained to lie on the line $ v_{1}=v_{2},$ and the linear components are constrained to be alternately horizontal and vertical.
Although the triangle above approximates the $w_{r}$ region, it also includes a region with $v_{2}>\chi _{\alpha }^{2}.$ The next theorem shows that no $CR_{LR}$ augmented with such a region, nor the CR suggested in VG2, can be size correct.
This might suggest a truncated version of this simple test with augmentation region
and critical region $\overline{CR}_{b}=CR_{LR}\cup w_{b}^{3}$, with superscript indicating that $w_{b}^{3}\ $ is contained in $A_{3}$. We show that there is an optimal $b,\bar{b},$ for this $\overline{LR}(b)$ test which is size-correct. We find $\bar{b}$ numerically such that NRP$<\alpha $ for $ 0\leq \lambda \leq \lambda _{0}$, using a relevant $\lambda _{0},$ and for $ \lambda >\lambda _{0}$ show that for this $\bar{b}$, the NRP is smaller than $\alpha $ using the following theorem.
Like the exact test, the $\overline{LR}(\bar{b})$ test is not coherent, however, as we show in Section (ref).
In view of the earlier comments, we consider the simple augmentation of $ CR_{LR}$ by the triangular region bounded from below by $v_{1}=bv_{2}$ and from above by $v_{1}=v_{2}<\chi _{\alpha }^{2}:$
The test is very simple: reject $H_{0}$ if $v_{1}\geq \chi _{\alpha }^{2}$ or $v_{1}/v_{2}>b$. It can be carried out by doing a LR\ test first. If it does not reject, then check if the ratio of the two test statistics is larger than the critical value. We denote the critical regions $CR_{LR}\cup w_{b}$ as $CR_{b}$ and refer to a test with given value of $b$ as an $LR(b)$ test. Obviously, $LR(1)$ is the LR test. The question to be addressed is how to choose the appropriate value of $b$. Some properties of this class of tests that can be used to choose $b$ are given next.
Define the discrepancy function as the difference between the NRP of the simply-augmented LR test and the nominal value $\alpha $:
Since the initial objective was to improve the behavior of the NRP near the origin, one possibility would be to choose the value of $b$ for which the test whose NRP is correct at the origin, i.e, the value satisfying $ D_{\alpha }(b,0)=0$ . This produces a test whose NRP is correct at $\lambda =0$, and also as $\lambda \rightarrow \infty ,$ but its NRP will be above the nominal level for intermediate values of $\lambda $.
The following result - partially reiterating Theorem (ref) - says that there is no member of this class of tests that has NRP equal to $\alpha $ (i.e., $D_{\alpha }(b,\lambda )=0)$ for all $\lambda $:
Now $D_{\alpha }(b,\lambda )$ is obviously continuous in $b$ on the interval $(0,1],$ and, for fixed $\lambda ,$ is (strictly) monotonic decreasing in $ b, $ with $D_{\alpha }(0,\lambda )=(1-2\alpha )G(\chi _{\alpha }^{2};\lambda )>0 $ when $\alpha <1/2,$ and $D_{\alpha }(1,\lambda )=-\alpha G(\chi _{\alpha }^{2};\lambda )<0$. This proves the following result:
Table (ref) illustrates the behavior of $ b(\lambda )$ for a few values of $\lambda $ when $\alpha =.05.\medskip $
\spacingset{1.2}
\spacingset{1.9} The NRP of $LR(b(\lambda ))$ does not exceed the nominal level $\alpha ,$ and choosing the largest $b(\lambda )$ results in NRP$<\alpha $. But this should hold for all values, also for\ $\lambda $ values not in Table (ref), even when $\lambda \rightarrow \infty $. The LR test has this asymptotic property, since $G(\chi _{\alpha }^{2};\lambda )\rightarrow 0$ as $\lambda \rightarrow \infty ,$ and this property is shared by all members of the class of $LR(b)$ tests:
So $\lim_{\lambda \rightarrow \infty }P_{CR_{b}}(\lambda )=\alpha $ for any $ b\in (0,1]$, including $P_{CR_{LR}}(\lambda )$, implying $P_{w_{b}}(\lambda )\rightarrow 0$ as $\lambda \rightarrow \infty $. Theorem (ref) implies nevertheless that there will exist $ \lambda $ such that $D_{\alpha }(b,\lambda )>0,$ because for any $b<1$, $ CR_{b}$ will include an area in $A_{2}$. The probability of this area will be very small however, if this $\lambda $ is large.
These considerations suggest that a reasonable approach to choosing $b$ would be to select the smallest value for which the maximum discrepancy as $ \lambda $ varies is bounded (small). This will produce a test with maximum power, subject to the constraint that overrejection is below the chosen bound.
Given that smaller $b$ implies higher power, but a value of $b$ too small leads to an invalid, over-sized test, we choose the smallest $b$ that is still (approximately) size correct. This is akin to choosing the largest $ b(\lambda )$ in Table (ref). So we choose the smallest $ b<1 $ such that $D_{\alpha }(b,\lambda )\leq \varepsilon ,$ where $ \varepsilon $ is a small number that needs to be chosen in the absence of a theoretical justification. It should be larger than 0 since $\varepsilon =0$ would imply correct size and $b=1,$ which is the LR test. Our choice $ \varepsilon =10^{-9}$ although theoretically positive, is negligible in practice and choosing $\varepsilon =10^{-16}$ only changes $b$'s, NRPs, and power by a small margin.
For each $\alpha $ percentile we determine the optimal (smallest) $b$ numerically such that $D_{\alpha }(b,\lambda )\leq \varepsilon $. Results are given in Appendix B and shown graphically in Figure (ref) for the truncated and the simply-augmented LR tests. For the simply-augmented LR test numerical values are given in Table (ref) at the end of the paper.
The optimal $b(\alpha )$ is monotonically decreasing in $\alpha $ for all three cases. For $\alpha <.05$ there is very little difference between the $ b $'s. For $\alpha \geq 1/2$ we have $b=0$ for the truncated version that then has CR $A_{1}\cup A_{3}$ and the NRP $=1/2\geq \alpha $ for all $ \lambda $.
Figure (ref) shows the $.05$ and $.10$ level NRPs for the LR and augmented LR tests. The VG2 test is not displayed, but would show a horizontal line, deviating from the level $\alpha $ by less than $10^{-9}$ for all $\lambda \geq 0$.
The LR starts at $\alpha ^{2}$ when $\lambda =0$ and increases to $\alpha $ as $\lambda \rightarrow \infty $. The NRP for $\overline{LR}(\bar{b})$ drops below that of $LR(b)$ for larger values of $\lambda $ when $w_{b}^{2}$ dominates $w_{\bar{b}}^{3}$ in probability terms. For small values of $ \lambda $, the NRP of the truncated version is larger, and progressively so with increasing $\alpha <1/2$. For $\alpha \geq 1/2$, the truncated version has $\bar{b}=0$ and is the exact similar test when $\alpha =1/2=NRP$ for all $\lambda $.
It is trivially true that the power of the $LR(b)$ and $\overline{LR}(\bar{b} )$ tests cannot be less than that of the LR test. The power function of the augmented LR test can be calculated in exactly the same way as we have done for the NRP. For the simply-augmented LR test and based on density ((ref)) this is:
the first term being the power of the LR test. The power function is obviously symmetric in $(\lambda _{1},\lambda _{2}).$ Again, as either noncentrality parameter goes to infinity, the power of the augmented test approaches that of the LR test. The power of the $LR(b)$ and $\overline{LR}( \bar{b})$ tests, together with that of the LR test itself (in brackets), are given in Table (ref) for a selection of values of $(\lambda _{1},\lambda _{2}).$ The table is, of course, symmetric.
\spacingset{1.2}
\spacingset{1.9} It is clear that, for $(\lambda _{1},\lambda _{2})$ near the origin, the LR test has poor power, and that the simply-augmented LR test of size $.05$ improves considerably upon it. The truncated version even more so since it has smaller $b$, so more area of $A_{3}$ near the origin is added. In Appendix C we display the power difference $P_{CR_{b}}(\lambda _{1},\lambda _{2})-P_{CR_{LR}}(\lambda _{1},\lambda _{2})$ for values of the $\lambda _{i}\in \lbrack 0,10].$ It is evident that the power difference is quite small for large $\lambda _{i},$ but substantial for $(\lambda _{1},\lambda _{2})$ near the origin.
We also show $P_{CR_{b}}(\lambda _{1},\lambda _{2})-P_{\overline{CR}_{\bar{b} }}(\lambda _{1},\lambda _{2})$ in Appendix B, which is much closer to zero and implies that the truncation restriction has little effect on the power. But the price being paid is lack of coherency which we will discuss next.
The LR critical region has the desirable property that if an observed point $ (v_{1},v_{2})$ falls in the rejection region at level $\alpha ,$ it also falls in the rejection region at every level larger than $\alpha .$ That is, the CR at a given level properly contains that at any smaller level. The LR test is thus coherent for inference on $H_{0}.$ An important property of the simply-augmented LR test as we have constructed it is that it retains this coherency property. This is perhaps best illustrated graphically. Figure (ref) shows the respective critical regions for the levels $\alpha =.01,.05,$ and $.1$: for the truncated version on the left, the simple version on the right. It is clear that the proposed approach provides coherent inference on $H_{0}$ in this sense, but only for the simply-augmented test, not the truncated version.
The coherence property just mentioned suggests that we can define a p-value for any observed point $(v_{1},v_{2})$ by reference to the critical regions $ CR_{b(\alpha )}.$ To do so, we simply determine the value of $\alpha ,$ say $ \alpha _{0},$ for which the observed point lies on the boundary of the critical region $CR_{b(\alpha _{0})}.$ Any point in the region $CR_{b(\alpha _{0})}$ lies on a boundary with a smaller level than $\alpha _{0},$ and in this sense are \textquotedblleft more extreme\textquotedblright\ under the null hypothesis than the observed point. The value $\alpha _{0}$ then has a natural interpretation as the p-value for the observed point.
To define $\alpha _{0}$ explicitly we make three observations: First, $ b(\alpha )$ is strictly monotonic and therefore has an inverse, so for each $ b\in \lbrack 0,1]$ there is a unique value $\alpha _{b}$ satisfying $ b=b(\alpha _{b})$. \ Second, each $v_{1}\geq 0$ yields a value $\alpha _{1}=1-G(v_{1})\in \lbrack 0,1].$ Third, every point $(v_{1},v_{2})\in V$ lies on either the horizontal part of the boundary of some $CR_{b(\alpha )},$ or on the sloping part. If it lies on the horizontal part, such that $ v_{1}/v_{2}<b(\alpha _{1}),$ then $\alpha _{0}=\alpha _{1}=1-G(v_{1})$. If it lies on the sloping part then $\alpha _{0}=b^{-1}(v_{1}/v_{2})$ and $ v_{1}/v_{2}=b(\alpha _{0}).$
Determining the p-value is now straightforward using Table (ref) : look up $v_{1}$ in the $\chi_{\alpha }^{2}$ column and note the corresponding $\alpha _{1}.$ If $v_{1}/v_{2}\leq b(\alpha _{1})$ then this $ \alpha _{1}$ is the p-value. If not, then $(v_{1},v_{2})$ must be on the sloping part of the boundary, so look up $v_{1}/v_{2}$ in the $b(\alpha )$ column. The p-value is the corresponding $\alpha $. One can interpolate for additional accuracy.
In order to demonstrate the applicability of the test and performance more broadly than the basic Gaussian setting, we carried out simulations using nonnormal and heteroskedastic disturbances, as well as a logit as a nonlinear model. Second, we include an empirical illustration from business.
The asymptotic normal distribution of two independent $t$-statistics forms the basis of the simply-augmented LR test. The approximation is valid for various estimation methods and error distributions asymptotically, but may be less accurate in small samples. Finite-sample critical values will deviate from those of the normal distribution, even in the case of homoskedastic normally distributed disturbances, when the two test statistics $f_{1}$ and $f_{2}$ are $F$-distributed with differing degrees of freedom, which further destroys the symmetry. To investigate the accuracy of the approximation in different circumstances, we used simulations to determine NRPs for different sample sizes and error distributions for $u_{1}$ and $u_{2}$: standard normal, student-$t$ (fat-tailed, with 5 degrees of freedom), chi-squared (skewed, with 3 degrees of freedom), and log-normal (skewed and fat-tailed), all standardized to mean zero and variance one. Sample sizes of 50, 100, 250, and 500 were considered, and $x_{i}\sim N(0,1)$ . Ordinary, rather than robust $t$-statistics are used. Given the NRP-deviation from $0.05$ in Figure (ref), the noncentrality parameter was chosen (approximately) as $\lambda _{2}\in \{0,1,...,25\}$ by solving $\lambda _{2}=\theta _{2}^{2}/(n-2)$ in (A.4) of Appendix A for $\theta _{2}$, and $\lambda _{1}=0$ implies $\theta _{1}=0$. Figure (ref) shows the NRP based on $P_{CB_{b}}(\lambda _{2})$, abbreviated to $P_{CB_{b}}$ in the figure's legend, and the simulated NRPs based on $10^{6}$ replications. For the normally distributed error case, the simulated NRPs are systematically higher than $P_{CB_{b}}$ with the largest deviation of approximately $0.006$ when $n=50$, in line with an $F$-distribution based critical value being larger than based on the $\chi ^{2}$ distribution. The simulated NRPs quickly converge to $P_{CBb},$ however, as the sample size increases with almost no deviance for $n=500$. In the case of errors that are standardized $t$-distributed or $\chi ^{2}$ -distributed, the results are very similar, although the convergence to $ P_{CB_{b}}$ is slower in the sample size. For the standardized log-normal error distribution, the maximum deviation is larger/smaller than for the other distributions for small/large values of the noncentrality parameter. Overall, the simulation results confirm that the approximation made in ((ref)) is very accurate, leading to a maximum overrejection of only $0.0025$ when $n=100$ for a wide variety of error distributions.
We continue to investigate the efficacy of employing robust $t$-statistics under homoskedasticity or in the presence of heteroskedasticity: $ (u_{1,i},u_{2,i})^{\prime }\sim N(0,\sigma (x_{i})I_{2})$ with (i) $\sigma (x_{i})=1$, (ii) $\sigma (x_{i})=|x_{i}|$, and (iii) $\sigma (x_{i})=\exp (0.4x_{i})$. Table (ref) shows that when ordinary $t$ -statistics are applied, the NRP may be as high as $20.8\%$ in case (ii). When robust $t$-statistics are used, the NRPs of the simply-augmented LR test are within $[4.6\%,6.5\%]$, while those of the LR\ test can be as low as $0.4\%$. The application of robust $t$-statistics only marginally increases the NRPs under homoskedasticity. Hence, in empirical research, it is recommended to use robust $t$-statistics. \spacingset{1.2}
\spacingset{1.9}
The simply-augmented LR test is applicable as long as the $t$-ratios are asymptotically standard normally distributed. This occurs more generally and one example of a nonlinear model with binary dependent variable is the following logit model with intercepts $\mu _{m}$ and $\mu _{y\ast }$:
where $\varepsilon _{i}\sim Logistic(0,1)$. The null of no mediation effect is still $H_{0}:\theta _{1}\theta _{2}=0$, which can be tested using the simply-augmented LR test, where $t_{1}$ and $t_{2}$ now denote the $t$ -ratios of $\hat{\theta}_{1}$ and $\hat{\theta}_{2}$ in the estimated regression ((ref)) and logit model ((ref)) respectively. The simulation closely follows the specification used in MacKinnon2007: $x$ is a dichotomous variable with equal numbers in each group, $u_{1}$ is either standard normal, student-$t$ (fat-tailed, with 5 degrees of freedom), or chi-squared (skewed, with 3 degrees of freedom) distributed, all standardized to mean zero and variance one. Parameter values for $\theta _{2}\in \{0.0,0.05,...,2.0\}$ were chosen to encompass the small (0.14), medium (0.39), large (0.59) and very large value (1.0) considered by MacKinnon2007, and $\theta _{1}=0$ (null hypothesis), $ \tau =0$ (irrelevant by invariance), and $\mu _{m}=\mu _{y\ast }=0$. The simulated NRPs in Figure (ref) are based on $10^{6}$ replications and are shown only for the normal and student-$t$ error distributions since the results for the chi-squared distribution are almost identical. The largest overrejection occurs in the normal case, with a maximal deviation of $0.006$ when $n=50$. For sample sizes $n\geq 100$, the new test is well-behaved with a maximum deviation of only $0.0026$. In the student-$t$ case, the deviations from the nominal level seem to be even smaller than in the normal case. For $n=500$ the figures are virtually identical to the results in Figure (ref) based on asymptotic theory. In summary, the simply-augmented $LR$ test also performs as expected in the important case of a binary dependent variable.
To illustrate our simply-augmented LR test, we use data from hayes2017introduction; this data set is called ESTRESS and can be downloaded from www.afhayes.com. The study involves entrepreneurs who were members of a networking group for small business owners; see pollack2012moderating. They answered an online survey about the recent performance of their business and also about their emotional and cognitive reactions to the economic climate. Hence, $Y$, $X$ and $M$ denote disengagement from entrepreneurial activities ($withdraw$), economic stress ( $estress$) and depressed affect ($affect$) respectively. There are also three confounding variables $C_{1}$, $C_{2}$ and $C_{3}$ that are related to entrepreneurial self-efficacy ($ese$), gender ($sex$) and length of time in the business ($tenure$); see Figure 4 of hayes2017introduction for a causal diagram. We focus on females with short $tenure$ (less than 0.6 years). OLS\ gives the following results (showing ordinary $t$-values in parentheses):\footnote{ Misspecification tests do not reject normality or homoskedasticity, with the following p-values:\newline Jarque-Bera (normality): 0.5078 ($M$) & 0.1848 ($Y$) and Breusch-Pagan (homoskedasticity): 0.2016 ($M$) & 0.1019 ($Y$).}
From these estimation results, we get the following test statistics: $ (f_{1},f_{2})=(1.277,1.254)$ and $(v_{1},v_{2})=(1.254,1.277)$, so that $ LR=1.254$. Since $LR<3.84$, the null of no mediation is not rejected at $.05$ level using the LR test or Sobel's test. Both $F$-statistics are very similar, however, leading to a ratio $v_{1}/v_{2}=.982$ that is larger than $ b(.05)=.8744$. The optimal simply-augmented test rejects, with a p-value of $ .00593$ using interpolation and Table (ref). The new test establishes a significant mediation effect that remains undetected by the conservative LR\ and Sobel tests.
We have demonstrated constructively that exact similar tests of the no-mediation hypothesis exist for tests of the nominal levels that are typically used in practice. Because the exact test described is incoherent, we have also proposed a basic modification of the LR test which is called the simply-augmented LR test. This new test to a large extent remedies the two main deficiencies of the standard LR test: very small NRP for small values of the nuisance parameter, and very poor power near the origin of the parameter space. This simply-augmented test is extremely easy to understand and apply in practice, and has the important property of coherence.
The earlier paper by VG2-2021, in which the idea of augmenting the $ LR$ critical region was first proposed for this problem, provides a more sophisticated augmentation method. It performs somewhat better in terms of NRP and power, but we show that it cannot be size correct.\ The test proposed here has the advantage of extreme simplicity and ease of implementation in empirical research. It requires only one additional number in conjunction with the standard critical value. Moreover, it is coherent, which is essential for deriving and interpreting p-values in applied work. We give results for each percentile in the table at the end of the paper, and show how p-values can be easily calculated and reported. A simulation study confirms the asymptotic approximations with non-Gaussian disturbances, with heteroskedastic disturbances, and in a relevant model for binary dependent variables.
There is no doubt that both the approach discussed in this paper, along with that of Van Garderen and Van Giersbergen, fall into the class of tests - departures from the likelihood ratio method - that is frowned upon by Perlman1999. Both lead to critical regions that imply rejection of the null hypothesis when the observed sample point is close to the origin. Moreover, the arguments claiming an improvement over the LR test are certainly based firmly on the Neyman-Pearson criteria of size and power. If one does not approve, the LR test is still available of course.
\spacingset{1.6}
\spacingset{1.2}
\numberwithin{equation}{section} \setcounter{section}{0} \setcounter{table}{0} \setcounter{figure}{0}
\spacingset{1.9}
{ \LARGEAppendix: Improved Tests for Mediation}
The testing problem for $H_{0}:\theta _{1}\theta _{2}=0,$ respects a number of symmetries. It allows restricting our attention for an optimal solution in terms of the maximal invariant statistic, see e.g. Davison2003, which we derive now under a Gaussianity assumption:
The statistics of interest are the sufficient statistics in the Gaussian model, which are simply the OLS estimates for ((ref)) and ((ref)), including the residual sums of squares, i.e. $(\hat{\theta} _{1},s_{11},\hat{\tau},\hat{\theta}_{2},s_{22})$ as defined in the proof of Theorem (ref). This theorem derives the maximal invariants which reduce the five-dimensional sufficient statistic to just two dimensions.
Proof. The maximal invariants $T_{1},T_{2}$ are the usual $t$ -statistics for testing $\theta _{1}=0$ and$\ \theta _{2}=0$. Under the Gaussianity assumption in equations ((ref)) and ((ref)), the MLEs for $\theta _{1}$ and $(\tau ,\theta _{2})$ are the OLS estimators and the MLEs for variances $\sigma _{11}$ and $\sigma _{22}$ are the residual sums of squares $s_{11}$ and $s_{22}$ divided by sample size $n$:
Here, for any matrix $A$ of full column rank, $M_{A}=I_{n}-A(A^{\prime }A)^{-1}A^{\prime }$. The distributions of the sufficient statistics are, respectively:
$s_{11}/\sigma _{11}\sim \chi ^{2}(n-1),$ $\hat{\theta}_{1}\sim N(\theta _{1},\sigma _{11}/s_{xx}),$ where $s_{xx}=x^{\prime }x,$ and $s_{22}/\sigma _{22}\sim \chi ^{2}(n-2)$. The joint density of the sufficient statistics under Gaussian assumptions may be written down directly from these facts, and is equivalent to the likelihood for $(\theta _{1},\sigma _{11},\tau ,\theta _{2},\sigma _{22}).$ The joint distribution of the sufficient statistics is a product of the form:
The transformations $s_{11}\mapsto a_{1}s_{11}$ and $s_{22}\mapsto a_{2}s_{22}$ with $a_{1},a_{2}>0$ leave the joint density of $\left( s_{11},s_{22}\right) $ in the same family with $\left( \sigma _{11},\sigma _{22}\right) $ replaced by $\left( a_{1}\sigma _{11},a_{2}\sigma _{22}\right) $ and have no bearing on the hypothesis under test. The same parameters $\sigma _{11}$ and $\sigma _{22}$ are present in the other components, so we need to transform the remaining variables accordingly, namely by: $\hat{\theta}_{1}\mapsto \sqrt{a_{1}}\hat{\theta}_{1}$ and $\hat{ \theta}_{2}\mapsto \sqrt{a_{2}/a_{1}}\hat{\theta}_{2}$. And, since $\tau $ is not involved in the inference problem, we may transform $\hat{\tau}$ by the affine transformation $\hat{\tau}\mapsto \sqrt{a_{2}}\left( \hat{\tau} +c\right) $.
These transformations preserve the family of distributions for the sufficient statistics (and MLEs), and the induced transformation on the mediation effect is that $\theta _{1}\theta _{2}\mapsto \sqrt{a_{2}}\theta _{1}\theta _{2}$. Thus, the transformations do not change the truth or falsity of the hypothesis under test (i.e. $H_{0}$ is true before iff it is true after the transformation). The transformations on $\hat{\tau}$ are transitive, so no invariant test can depend on $\hat{\tau}.$ We can therefore restrict attention to the four remaining statistics\ $ \left( \hat{\theta}_{1},s_{11},\hat{\theta}_{2},s_{22}\right) $, and the group $\mathbf{K}$, say, of (scale) transformations of them. The invariance of $(T_{1},T_{2})$ under the transformations is obvious. To show that $ (T_{1},T_{2})$ are maximal we need to show that $T_{1}(\hat{\theta}_{1},\hat{ \theta}_{2},s_{11},s_{22})=T_{1}(\tilde{\theta}_{1},\tilde{\theta}_{2}, \tilde{s}_{11},\tilde{s}_{22})$ and $T_{2}(\hat{\theta}_{1},\hat{\theta} _{2},s_{11},s_{22})=T_{2}(\tilde{\theta}_{1},\tilde{\theta}_{2},\tilde{s} _{11},\tilde{s}_{22})$ implies that there exists a group element $K\in \mathbf{K}$ such that $(\tilde{\theta}_{1},\tilde{\theta}_{2},\tilde{s}_{11}, \tilde{s}_{22})=K(\hat{\theta}_{1},\hat{\theta}_{2},s_{11},s_{22})$.
Thus, assume that
and
Then $\tilde{\theta}_{1}=\sqrt{a_{1}}\hat{\theta}_{1}$ with $a_{1}=\tilde{s} _{11}/s_{11}$, and $\tilde{\theta}_{2}=\sqrt{a_{2}/a_{1}}\hat{\theta}_{2}$ with $a_{2}=\tilde{s}_{22}/s_{22}$. Since also $\tilde{s}_{11}=a_{1}s_{11}$, and $\tilde{s}_{22}=a_{2}s_{22}$, this shows that the invariance of $ (T_{1},T_{2})$ implies that the two sets of statistics are related by a group element, so $(T_{1},T_{2})$ are indeed maximal. The same argument applies to the induced group acting on the parameter space, and the last statement is a well-known property of maximal invariants.
Proposition 1 can be proved using the assumptions and theorems in White1980 as is explicitly done in VG2-2021.
The noncentral $\chi _{\kappa }^{2}$ density, $g_{\kappa }(f;\lambda ),$ with noncentrality $\lambda $, plays a central role throughout this paper and its appendix.\ It can be expressed in several ways, including the Poisson mixture exploited in the proof of Proposition 5:
where $g_{\kappa }(f)=[2^{\frac{\kappa }{2}}\Gamma (\frac{\kappa }{2} )]^{-1}\exp \{-\frac{1}{2}f\}f^{\frac{\kappa }{2}-1}$ denotes the $\chi _{\kappa }^{2}$ density function, and we write $g_{1}(f)$ simply as $g(f).$ The corresponding CDFs are denoted by $G_{\kappa }(\cdot ;\lambda ),$ and $ G_{\kappa }(\cdot )$ in the central case, the subscript being omitted when $ \kappa =1.$ For $\alpha \in \lbrack 0,1],$ we define $\chi _{\alpha }^{2}$ by $G(\chi _{\alpha }^{2})=1-\alpha .$
From Equation (6) in Vaughan1972, the joint density of the order statistics for $(v_{1},v_{2})\in V=\{(v_{1},v_{2});0\leq v_{1}\leq v_{2}<\infty \}$ is as given in part (i) of the following proposition, which also gives complete details of the distribution of the order statistics$: \footnote{ Since the order statistics are maximal invariants under the action of the symmetric group $S_{2}$ on $(f_{1},f_{2}),$ the main result in part (i) can also be obtained by invoking Stein's method of obtaining the density of the maximal invariant by averaging the joint density over the group.}$
Proof of Proposition (ref) and remarks
Part (i) is direct from Vaughan1972. Part (iii) is simply the fact that, on integrating over $v_{1}<v_{2}<\infty ,$ we have
Remarks
(i) We can specialize these results for the null case when one noncentrality parameter vanishes
where $\lambda $ is the nonzero noncentrality parameter.
(ii) It is trivial to check that the derivatives of $H(v;\lambda )$ and $ H(v_{1};\lambda _{1},\lambda _{2})$ yield the densities given in ((ref)) and ((ref)).
(iii) Observe that
This and similar identities are useful, for example, to verify that the joint density integrates to one:
It is not difficult to obtain the following probabilities under the null: for any $z>0$,
Excluding the region $A_{r}(z)$ from $A_{3}$ leaves $r+1$ disjoint triangles lying along the $45^{\circ }$ line, the region $w_{r}(z)\subset A_{3},$ and it is easy to see that the null probability content of this region is
This differs from the target value $\alpha G(\chi _{\alpha }^{2};\lambda )$ by
This is a linear combination of $r+1$ noncentral chi-square CDFs, the $ G(z_{i};\lambda )$ for $i\in \{1,...,r\}$, and $G(z_{r+1};\lambda )=G(\chi _{\alpha }^{2};\lambda ),$ and vanishes for all $\lambda $ if and only if the coefficients of all $r+1$ terms that involve $\lambda $ vanish. With $ \alpha =1-G(z_{r+1}),$ these conditions give rise to a system of $r+1$ linear equations in $r+1$ unknowns $G(z_{i}),i=1,...,r+1$ (we take $z_{0}=0,$ so that $G(z_{0})=0),$ and these determine $z_{1},...,z_{r}$ and $ z_{r+1}=\chi _{\alpha }^{2},$ hence $\alpha $. \ It is easy to see that the matrix of the system is non-singular, so the solution is unique, and it is straightforward to check that the solution is:
so that $\alpha =1-G(z_{r+1})=(r+2)^{-1}.$
Define a rectangle $R$ inside $AR_{LR}$ by first defining $\chi _{\alpha +\epsilon }^{2}=z_{1}$ and $\chi _{\alpha -\epsilon }^{2}=z_{2}$ for $ 0<\epsilon \leq \alpha \leq \frac{1}{2},$ such that $1-G\left( z_{1}\right) =\alpha +\epsilon ,$ $1-G\left( z_{2}\right) =\alpha -\epsilon $, $G(\chi _{\alpha }^{2})-G(z_{1})=\epsilon ,$ $G(z_{2})-G(\chi _{\alpha }^{2})=\epsilon ,$ $z_{1}\leq \chi _{\alpha }^{2}\leq z_{2},$ and then
in the top left-hand corner of $A_{2}$. The probability content of $R$ is:
The rejection probability of the LR\ test augmented with $R$ is:
Hence the test is correctly sized iff
or
We now prove that this inequality cannot hold by constructing a lower-bound and showing that it exceeds $\frac{\alpha }{\epsilon }$ for any $\alpha $ and $\epsilon $, as $\lambda \rightarrow \infty .$ The bound is constructed by deriving a lower bound for the numerator and an upper bound for the denominator. First note that the CDF\ of the noncentral Chi-square distribution with one degree of freedom can be written as\footnote{ By a simple transformation of the random variable $Z=X^{2}$ with $X\sim N( \sqrt{\lambda },1)$.}
where $\Phi (\cdot )$ denotes the cdf of the standard normal $N(0,1)$ and the second line is written such that the argument of $\Phi (\cdot )$ is positive when $\lambda >z$. The numerator in ((ref)) can therefore be written as
Using the Hermite-Hadamard inequality for convex $f\left( x\right) ,$ which states:
and given that the pdf of the standard normal $\phi (x)$ is convex for $x>1$ we obtain:
So for $z_{2}>\chi _{\alpha }^{2}>z_{1}>1$ we can write for the four terms in ((ref))
since $\phi (\sqrt{\lambda }+1/2\sqrt{z_{1}}+1/2\sqrt{z_{2}})>0$. When $ \epsilon $ is small, $z_{1}$ and $z_{2}$ are very similar and the lower bound based on ((ref)) will be tight. \newline For the upper bound of the denominator $G(\chi _{\alpha }^{2};\lambda )$ we use the following inequalities:
see, e.g. Aggarwal2019 in Equation ((ref)):
Combining the bounds ((ref)) and ((ref) ) gives:
The scaling factor in front of the exponential function is positive and tends to infinity as $\lambda \rightarrow \infty $ since $z_{2}>\chi _{\alpha }^{2}>z_{1}$. The limit behavior of the exponential function is determined by the coefficient of $\sqrt{\lambda }$, which is positive (and real) since
where $\sqrt{z_{1}},\sqrt{z_{2}}$ and $\sqrt{\chi _{\alpha }^{2}}$ are the quantiles of a standard normal distribution for $1-(\alpha +\epsilon )/2,$ $ 1-(\alpha -\epsilon )/2,$ $1-\alpha /2$ respectively and the quantile function of a standard normal is convex for probabilities larger than 1/2. This implies that the ratio in ((ref)) can be made arbitrarily large by choosing $\lambda $ sufficiently large. This violates the condition that it should be smaller than $\alpha /\epsilon .$ Hence there is a $ \lambda _{0}$ such that $NRP(\lambda )>\alpha $ for all $\lambda _{0}<\lambda <\infty $ and the test is over-sized for finite $\lambda $. Of course, $\lim_{\lambda \rightarrow \infty }NRP(\lambda )=\alpha $ as shown in Proposition 7 of the paper.
We now consider the truncated-simply-augmented LR test defined by augmentation of $CR_{LR}$ with the region $w_{b}^{3}$. This region is defined as:
for $0<b\leq 1$.
The probability content of region $w_{b}^{3}$ is given by:
Note the change in the order of integration in the third line. Since the probability of $CR_{LR}$ is given by $\alpha -\alpha G(\chi _{\alpha }^{2},\lambda )$, the difference between level $\alpha $ and the NRP is given by $-\alpha G(\chi _{\alpha }^{2},\lambda )$. Noting that $G(\chi _{\alpha }^{2},\lambda )=\int_{0}^{b\chi _{\alpha }^{2}}g(v;\lambda )dv+\int_{b\chi _{\alpha }^{2}}^{\chi _{\alpha }^{2}}g(v;\lambda )dv$, the difference between level $\alpha $ and the NRP\ of the truncated-augmented critical region is given by the following discrepancy function:
where the two functions $f_{1}(\cdot )$ and $f_{2}(\cdot )$ are independent from $\lambda $. For a given $\lambda $, $P[w_{b}^{3}]$ increases as $b$ decreases. Hence, for each $\alpha $, there is an optimal value $b^{\ast }$ such that $D_{\alpha }(b^{\ast },\lambda )$ is numerically close to 0 for some value of $\lambda $. Using the result
of Cohen1988, the derivative of $\bar{D}_{\alpha }(b,\lambda )$ with respect to $\lambda $ is given by:
For $\alpha =.05$, an iterative procedure between choosing $b$ and finding the maximum value of $\bar{D}_{.05}(b,\lambda )$ gives $b^{\ast }=0.86978984806$, where the maximum is obtained at $\lambda ^{\ast }=3.844989224948$ such that $\bar{D}_{.05}(b^{\ast },\lambda ^{\ast })=-2.95336772\cdot 10^{-13}$. The figure below shows $\bar{D}_{.05}(b^{\ast },\lambda )$ for $1\leq \lambda \leq 25$ (left) and a close-up for $\lambda $ around $\lambda ^{\ast }$ (right).
If $\bar{D}_{.05}(b^{\ast },\lambda )\leq 0$ for all $\lambda \geq 0$, then the test based on $CR_{LR}$ augmented with region $w_{b}^{3}$ is a valid $.05 $-level test. The figure above numerically shows that $\bar{D}_{.05}(b^{\ast },\lambda )<0$ for $\lambda \leq 25$, but this does not guarantee that $\bar{ D}_{.05}(b^{\ast },\lambda )$ remains below zero for larger values of $ \lambda $. In Theorem 3 of the paper, it is claimed that $\bar{D}_{\alpha }(b^{\ast },\lambda )\leq 0$, so $\bar{D}_{\alpha }(b^{\ast },\lambda )$ is bounded from above by a function that becomes negative for sufficiently large $\lambda (\alpha )$. This bounding function is based on an upper bound of the two integrals in ((ref)) in such a way that relatively simple expressions in the pdf of the standard normal distribution $\phi (\cdot )$ are obtained. In the figure below, $f_{1}(v;b^{\ast })$ and $ f_{2}(v;b^{\ast })$ are shown in blue and orange respectively for $\alpha =.05$.
We are now in the position to prove the following claim that implies Theorem 3 of the paper:
For $\alpha <0.5$, there is a bounding function $ub_{\alpha }(\lambda )\leq 0$\ and a value $\lambda _{0}(\alpha )$\ such that $D_{\alpha }(b^{\ast },\lambda )\leq ub_{\alpha }(\lambda )$\ for $\lambda _{0}(\alpha )<\lambda <\infty $.
To show that the sum of the integrals in ((ref)) is negative, $f_{1}(v;b^{\ast })$ and $f_{2}(v;b^{\ast })$ are bounded from above. Note that for large values of $\lambda $, the pdf of the noncentral chi-square behaves very similar to the pdf of a scaled standard normal, i.e.
which is an exponentially increasing function in (the scalar) $v$. As $ \lambda $ increases, more weight is given to larger $v$-values. For $\alpha =.05$, the bounding function is shown as the dashed red line in the figure and corresponds to the function $\gamma _{0}+\gamma _{1}\sqrt{v}$ with $ (\gamma _{0},\gamma _{1})=(0.07823,-0.04917)$. Note that the bound is only tight for larger values of $v$, but this is sufficient as this area gets the highest weight. The bounding function is chosen because the antiderivative of $\sqrt{v}g(v;\lambda )$ is relatively simple:
Note that in the first line, the arguments of $G(\cdot )$ are switched, which is done intentionally and is not a typo. Hence, the integrals in ((ref)) are bounded by:
Substitution of ((ref)) in ((ref)) leads to
where differences $\Phi (b)-\Phi (a)$ are such that $b>a$. In order to investigate the behavior of $bound(\lambda )$ as $\lambda $ becomes large, each of the three terms within parenthesis in ((ref))-((ref)) will be bounded. For the bounding function, it is important to note that $ \gamma _{0}>0$, whereas $\gamma _{1}<0$.
For the first two terms, the inequalities shown in ((ref)) are used. Using the upper bound for $\Phi (\sqrt{\lambda }+\sqrt{\chi_{\alpha }^{2}})$ and the lower bound for $\Phi (\sqrt{\lambda }-\sqrt{\chi_{\alpha }^{2}})$, the first term in ((ref)) is bounded by
Taking appropriate upper and lower bounds, the second term is bounded by
The last term is the smallest in magnitude since it does not involve integrals over $\phi $. For $\lambda >\chi_{\alpha }^{2}/4$, we have $\phi ( \sqrt{\chi_{\alpha }^{2}}-\sqrt{\lambda })>\phi (\sqrt{\lambda })>\phi ( \sqrt{\chi_{\alpha }^{2}}+\sqrt{\lambda })$, so the first term in ((ref)) can be bounded by:
Using the three upper bounds, the total can be bounded from above by
Algebraic simplifications show that the upper bound $ub(\lambda )$ can be written as
where
When $\gamma _{0}>0$ and $\gamma _{1}<0$, $c_{2,\alpha }(\lambda )>0$ for $ \lambda >0$ and $\exp (-\chi _{\alpha }^{2}/2-\sqrt{\chi _{\alpha }^{2}\lambda })>0$ tends exponentially to zero as $\lambda \rightarrow \infty $. Hence, we see that the sign of $ub_{\alpha }(\lambda )$ is mainly determined by the sum $\exp (-\chi _{\alpha }^{2}/2+\sqrt{\chi _{\alpha }^{2}\lambda })c_{1,\alpha }(\lambda )-2\gamma _{1}$. Since $\exp (-\chi _{\alpha }^{2}/2+\sqrt{\chi _{\alpha }^{2}\lambda })>0$ becomes exponentially large as $\lambda $ increases, the sign of the sum will be determined by the sign of $c_{1,\alpha }(\lambda )$. The coefficient $ c_{1,\alpha }(\lambda )$ can be written as a ratio of two polynomials in $ \sqrt{\lambda }$, which is done in ((ref)). Since
the function $c_{1,\alpha }(\lambda )$ will become and stay negative for sufficiently large $\lambda $ when $\gamma _{0}<-\sqrt{\chi _{\alpha }^{2}} \gamma _{1}$. For selected values of $\alpha $, Table (ref) shows the optimal values $b^{\ast }$, the value for which the discrepancy function has a local maximum $\lambda ^{\ast }$ and the numerical value of $(\gamma _{0},\gamma _{1})$ such that $(\gamma _{0}+\gamma _{1}\sqrt{v})$ is larger than $f_{1}(v;b^{\ast })$ for $0\leq v\leq b^{\ast }\chi _{\alpha }^{2}$ and $f_{2}(v;b^{\ast })$ for $b^{\ast }\chi _{\alpha }^{2}\leq v\leq \chi _{\alpha }^{2}$. Furthermore, the last column shows $\lambda _{0}(\alpha )$ such that $ub_{\alpha }(\lambda )$ is below zero for $\lambda >\lambda _{0}(\alpha )$. For $\alpha =.05$, the ratio in ((ref)) given by
which turns negative for $\lambda >33.580,$ whereas Table (ref) shows that $ub_{.05}(\lambda )<0$ for $\lambda >33.64$ . Both functions converge to zero from below as $\lambda \rightarrow \infty . $ The figure below shows $\bar{D}_{.05}(b^{\ast },\lambda )$ as well as the upper bound $ub_{.05}(\lambda )$ in the region of $\lambda $ around $ \lambda _{0}(.05)$.
\spacingset{1.0}
\spacingset{1.9}
Under $H_{0}$ the probability content of the augmenting region is given by
Substituting for the density and evaluating the integral over $v_{1}$ produces, after simplification,
Integrating each term in the second line by parts gives
Then, transforming to $v=bv_{2}$ in the integral, we obtain
and the result follows.
Expanding the two noncentral components in the integrand in $A_{\alpha }(b;\lambda )$ as Poisson mixtures we have
and also
The two power series coincide for all $\lambda $ if\ and only if all coefficients agree, that is
for all $j.$ There is no $b\in (0,1]$ satisfying this equation for all $j.$
It is well-known that $\lim_{\lambda \rightarrow \infty }G(z;\lambda )=0$ for any finite $z>0.$ The term $A_{\alpha }(b;\lambda )$ in $D_{\alpha }(b,\lambda )$ is evidently positive, and is less than $G(\chi _{\alpha }^{2}/b)G(\chi _{\alpha }^{2};\lambda )+G(\chi _{\alpha }^{2})G(\chi _{\alpha }^{2}/b;\lambda )$ for all $\lambda .$ Since this converges to $0$ as $\lambda \rightarrow \infty $ for any $b>0,$ both terms in $D_{\alpha }(b,\lambda )$ go to zero as $\lambda \rightarrow \infty .$
Below, an alternative expression is derived for the discrepancy function for the simply-augmented $LR(b)$ test with the same structure as in ((ref)). The probability content of region $w_{b}=w_{b}^{2}\cup w_{b}^{3}$ is given by
The discrepancy function is given by $D_{\alpha }(b,\lambda )=P_{w_{b}}(\lambda )+P_{A_{1}}(\lambda )-\alpha =P_{w_{b}}(\lambda )-\alpha G(\chi _{\alpha }^{2};\lambda )$. Noting that $-\alpha G(\chi _{\alpha }^{2};\lambda )=\int_{0}^{\chi _{\alpha }^{2}}-\alpha g(v;\lambda )dv$, we obtain
where the two functions $h_{1}(\cdot )$ and $h_{2}(\cdot )$ are independent from $\lambda $. Using the result in ((ref)), the derivative of $ D_{\alpha }(b,\lambda )$ with respect to $\lambda $ is given by:
When $b=1$, the area $w_{b}$ reduces to zero and the discrepancy function is negative for any $\lambda \geq 0$, i.e. $D_{\alpha }(1,\lambda )=P_{A_{1}}(\lambda )-\alpha =-\alpha G(\chi _{\alpha }^{2};\lambda )<0$. For $b<1$, the discrepancy function can become positive. The goal is to determine $b$ as small as possible given a maximum positive deviation. The shape of the discrepancy function changes with $\alpha $ and \thinspace $b$. For instance, for $\alpha =.05$, $D_{.05}(b,\lambda )$ has three stationary points for $0.86635\leq b\leq 0.87547$, but only one stationary point outside this interval given a grid of $b$-values. However, for $\alpha =.01$ , $D_{.01}^{\prime }(b,\lambda )$ has one root and hence $D_{.01}(b,\lambda ) $ only possesses one maximum for $b<1$; see Figure (ref) for $\alpha =.05$.
To find the smallest value of $b$ such that $D_{\alpha }(b,\lambda )<\epsilon $ for $\epsilon \in \{10^{-9},10^{-16}\}$ the following algorithm is used. For a given value of $b$, the number of sign changes in the derivative of $D_{\alpha }(b,\lambda )$ is determined for a grid of $\lambda $-values, i.e. $\lambda \in \{0.0001:0.01:5\}\cup \{5.2:0.2:30\}\cup \{31:1:150\}$ where ${\{a:b:c\}}$ is a regular grid between $a$ and $c$ with grid spacing $b$. For each sign change, a root of $D_{\alpha }^{\prime }(b,\lambda )$ is found by bisection using formula ((ref) ). The number of roots varies with $b$, and we identify two cases: case 1: only one root $\lambda _{1}^{\ast }$ is found or case 2: three roots $ \lambda _{1}^{\ast }\geq \lambda _{2}^{\ast }\geq \lambda _{3}^{\ast }$ are found. In each step, $b$ is reduced by $\Delta =0.0001$, i.e. $b\leftarrow b-\Delta $. This is stopped when, in case 1: $D_{\alpha }(b,\lambda _{1}^{\ast })>\epsilon $ or, in case 2: $D_{\alpha }(b,\lambda _{1}^{\ast })>\epsilon $ or $D_{\alpha }(b,\lambda _{3}^{\ast })>-\epsilon $. Next, $ \Delta $ is divided by 10, i.e. $\Delta \rightarrow \Delta /10$ and it is checked if for the slightly higher value of $b$, i.e. $b+\Delta $, in case 1: if $D_{\alpha }(b+\Delta ,\lambda _{1}^{\ast })<\epsilon $ or in case 2: if $D_{\alpha }(b+\Delta ,\lambda _{1}^{\ast })<\epsilon $ and $D_{\alpha }(b+\Delta ,\lambda _{3}^{\ast })<-\epsilon $. When this is not the case, repeat the procedure with an appropriate starting value, i.e. $b+10\Delta $ and the reduced $\Delta $ obtained before. All computations were done in Julia version 1.8.5, see Bezanson2017, and verified in Mathematica 13.1, see Mathematica. Table (ref) shows the results for $\alpha \in \{.01,.05,.1\}$, where for $\alpha =.01$ the largest root is determined outside the mentioned grid of $\lambda $ values using arbitrary precision arithmetic in Julia using setprecision(96) corresponding to approximately 32 significant digits.
\spacingset{1.0}
\spacingset{1.9}
Note the difference in scale between Figure (ref) and Figure (ref).
\spacingset{1.0}
\spacingset{1.0}