Tests for qualitative features in the random coefficients model

TL;DR

This paper develops statistical tests to detect qualitative features like modes and density increases in the joint distribution of random coefficients in a linear model, addressing the challenges of ill-posedness and high dimensionality.

Contribution

It introduces a multiple testing framework for identifying shape features of the joint density in random coefficients models, using Gaussian approximation techniques.

Findings

Effective detection of modes and density increases in simulated data

Method applied successfully to real data examples

Provides confidence statements for qualitative features

Abstract

The random coefficients model is an extension of the linear regression model that allows for unobserved heterogeneity in the population by modeling the regression coefficients as random variables. Given data from this model, the statistical challenge is to recover information about the joint density of the random coefficients which is a multivariate and ill-posed problem. Because of the curse of dimensionality and the ill-posedness, pointwise nonparametric estimation of the joint density is difficult and suffers from slow convergence rates. Larger features, such as an increase of the density along some direction or a well-accentuated mode can, however, be much easier detected from data by means of statistical tests. In this article, we follow this strategy and construct tests and confidence statements for qualitative features of the joint density, such as increases, decreases and modes. We propose a multiple testing approach based on aggregating single tests which are designed to extract shape information on fixed scales and directions. Using recent tools for Gaussian approximations of multivariate empirical processes, we derive expressions for the critical value. We apply our method to simulated and real data.

Tables

Table 1. Table 1: Simulated level and power of the test ( 4.1 ) for a mode considering uniform design 𝐗 i ∼ Unif [ − 5 , 5 ] 3 \mathbf{X}_{i}\sim\operatorname{Unif}[-5,5]^{3} and normal design 𝐗 i ∼ 𝒩 ​ ( ( 3 , 0 , 0 ) ⊤ , 2 ​ I ) similar-to subscript 𝐗 𝑖 𝒩 superscript 3 0 0 top 2 𝐼 \mathbf{X}_{i}\sim\mathcal{N}((3,0,0)^{\top},2I) . Results with theoretical quantiles are in the second column. Results where κ n 𝐭 j , h 0 , 𝐯 j ​ ( α ) superscript subscript 𝜅 𝑛 subscript 𝐭 𝑗 subscript ℎ 0 subscript 𝐯 𝑗 𝛼 \kappa_{n}^{\mathbf{t}_{j},h_{0},\mathbf{v}_{j}}(\alpha) in ( 4.1 ) are replaced by calibrated quantiles are in all other columns.

	${𝐗}_i\sim \mathrm{Unif}{[-5,5]}^3$			${𝐗}_i\sim 𝒩({(3,0,0)}^{\top },2I)$
$n$	power	level (cal.)	power (cal.)	level (cal.)	power (cal.)
250	9	4.5	91.0	4.7	66.1
500	15.7	4.6	99.1	4.5	78.8
1000	79.1	5.0	100	4.5	90.6

Table 2. Table 2: Same as Table 1 but for random coefficients model with intercept and ( X i , 2 , X i , 3 ) ⊤ superscript subscript 𝑋 𝑖 2 subscript 𝑋 𝑖 3 top (X_{i,2},X_{i,3})^{\top} from a standard bivariate Cauchy distribution. In the fourth and fifth column we assume that the density f 𝚯 subscript 𝑓 𝚯 f_{\bm{\Theta}} is known.

	$f_{𝚯}$ unknown		$f_{𝚯}$ known
$n$	level (cal.)	power (cal.)	level (cal.)	power (cal.)
250	4.8	91.3	5.0	93.3
500	5.2	99.0	5.2	99.7
1000	5.3	100	5.3	100

Table 3. Table 3: Same as Table 1 but for random coefficients model with intercept and ( X i , 2 , X i , 3 ) ⊤ ∼ 𝒩 ​ ( ( 0 , 0 ) ⊤ , I ) similar-to superscript subscript 𝑋 𝑖 2 subscript 𝑋 𝑖 3 top 𝒩 superscript 0 0 top 𝐼 (X_{i,2},X_{i,3})^{\top}\sim\mathcal{N}((0,0)^{\top},I) (second and third column) and ( X i , 2 , X i , 3 ) ⊤ ∼ Unif [ − 5 , 5 ] 2 (X_{i,2},X_{i,3})^{\top}\sim\operatorname{Unif}[-5,5]^{2} (fourth and fifth column).

	${(X_{i,2},X_{i,3})}^{\top }\sim 𝒩(0,1)$		${(X_{i,2},X_{i,3})}^{\top }\sim \mathrm{Unif}{[-5,5]}^2$
$n$	level (cal.)	power (cal.)	level (cal.)	power (cal.)
250	5.1	88.3	5.1	64.8
500	5.5	98.3	5.1	78.1
1000	5.1	100	5.5	89.4

Table 4. Table 4: Power of the multiscale test for modality. The mode tests for the bandwidths h = 2.5 ℎ 2.5 h=2.5 , h = 1 ℎ 1 h=1 and h = 0.5 ℎ 0.5 h=0.5 are in the second, third and fourth column. The combined mode test using h = 1 ℎ 1 h=1 and h = 0.5 ℎ 0.5 h=0.5 is in the fifth column.

$n$	$h=2.5$	$h=1$	$h=0.5$	$h\in \{0.5,1\}$
2000	100	0	0	25
5000	100	0	1	86.3
15000	100	0	68.7	100

Table 5. Table 5: Results of OLS for model Y = 𝛄 ⊤ ​ 𝐗 + 𝐮 ⊤ ​ 𝐗 𝑌 superscript 𝛄 top 𝐗 superscript 𝐮 top 𝐗 Y={\bm{\gamma}}^{\top}\mathbf{X}+\mathbf{u}^{\top}\mathbf{X} for d = 3 𝑑 3 d=3 and ( X i , 2 , X i , 3 ) ⊤ superscript subscript 𝑋 𝑖 2 subscript 𝑋 𝑖 3 top (X_{i,2},X_{i,3})^{\top} drawn from a standard bivariate Cauchy distribution. The distribution of the true random coefficients 𝛃 i subscript 𝛃 𝑖 {\bm{\beta}}_{i} is a standard Gaussian in column (i), the normal mixture 0.5 ​ 𝒩 ​ ( 𝟎 , 0.1 ​ I ) + 0.5 ​ 𝒩 ​ ( ( 2 , 0 , 0 ) ⊤ , 0.1 ​ I ) 0.5 𝒩 0 0.1 𝐼 0.5 𝒩 superscript 2 0 0 top 0.1 𝐼 0.5\mathcal{N}(\bm{0},0.1I)+0.5\mathcal{N}((2,0,0)^{\top},0.1I) in column (ii) and an exponential-2 distribution for β i , 1 subscript 𝛽 𝑖 1 \beta_{i,1} and ( β i , 2 , β i , 3 ) ⊤ ∼ 𝒩 ​ ( 𝟎 , 0.1 ​ I ) similar-to superscript subscript 𝛽 𝑖 2 subscript 𝛽 𝑖 3 top 𝒩 0 0.1 𝐼 (\beta_{i,2},\beta_{i,3})^{\top}\sim\mathcal{N}(\bm{0},0.1I) independent of the first component in column (iii).

	(i) ($n=1000$)		(ii) ($n=1000$)		(iii) ($n=1000$)
${\hat{\gamma }}_1$	0.00	(0.02)	1.00	(0.02)	2.03	(0.03)
${\hat{\gamma }}_2$	0.02	(0.02)	0.00	(0.01)	-0.02	(0.02)
${\hat{\gamma }}_3$	0.00	(0.02)	0.01	(0.01)	0.02	(0.02)

Table 6. Table 6: Results of OLS for model ( 4.2 ) with fixed coefficients.

${\tilde{\beta }}_1$	${\tilde{\beta }}_2$	${\tilde{\beta }}_3$
$0.1940$	$-0.0587$	$0.0005$
$\left(0.000\right)$	$\left(0.001\right)$	$\left(0.001\right)$

Equations

$$Y_i={\beta }_{i,1}{X}_{i,1}+{\beta }_{i,2}{X}_{i,2}+\dots +{\beta }_{i,d}{X}_{i,d},i=1,\dots ,n.$$

$$Y_i={\beta }_{i,1}+{\beta }_{i,2}{X}_{i,2}+\dots +{\beta }_{i,d}{X}_{i,d},i=1,\dots ,n,$$

$$Y_i=Q_i\cos ({\Phi }_i)+P_i\sin ({\Phi }_i),i=1,\dots ,n,$$

$$Rf(s,\mathbf{\theta })={\int }_{⟨\mathbf{b},\mathbf{\theta }⟩=s}f(\mathbf{b})d{\mu }_{d-1}(\mathbf{b})$$

$$S_i:=\frac{Y_i}{\Vert {\mathbf{X}}_i\Vert },{\mathbf{\Theta }}_i:=\frac{{\mathbf{X}}_i}{\Vert {\mathbf{X}}_i\Vert },i=1,\dots ,n.$$

$$F_{S|\mathbf{\Theta }}(x|\mathbf{\theta })=\mathbb{P}(S_1\le x|{\mathbf{\Theta }}_1=\mathbf{\theta })=\mathbb{P}(⟨{\mathbf{\Theta }}_1,{\mathbf{\beta }}_1⟩\le x|{\mathbf{\Theta }}_1=\mathbf{\theta })={\int }_{-\infty }^x(R{f}_{\mathbf{\beta }})(s,\mathbf{\theta })ds,$$

$$f_{S|\mathbf{\Theta }}(s|\mathbf{\theta })=(R{f}_{\mathbf{\beta }})(s,\mathbf{\theta }).$$

$$\Lambda f(s,\mathbf{\theta })={\mathcal{H}}_d{\partial }_s^{d-1}Rf(s,\mathbf{\theta }),$$

$$\mathcal{H}_{d}f(u)=\frac{1}{\pi}\lim_{\epsilon\rightarrow 0^{+}}\int_{(-\infty,u-\epsilon]\cup[u+\epsilon,\infty)}\frac{f(s)}{u-s}ds=\frac{1}{\pi}\text{p.v.}\int_{-\infty}^{\infty}\frac{f(s)}{u-s}ds\quad\big{(}u\in\mathbb{R}\big{)}$$

$$\phi =c_d^{-1}{R}^{\ast }\Lambda \phi ,$$

$${\int }_{{\mathbb{R}}^d}{\varphi }_{\mathbf{t},h}(\mathbf{b}){\partial }_{\mathbf{v}}{f}_{\mathbf{\beta }}(\mathbf{b})d\mathbf{b}=\sum _{k=1}^d{\mathrm{v}}_k{\int }_{{\mathbb{R}}^d}{\varphi }_{\mathbf{t},h}(\mathbf{b})\frac{\partial }{\partial {\mathrm{b}}_k}{f}_{\mathbf{\beta }}(\mathbf{b})d\mathbf{b}$$

$$\displaystyle\phi_{\mathbf{t},h}(\cdot)=\frac{1}{h^{d}\operatorname{Vol}(\mathbb{S}^{d-2})}\phi\Big{(}\frac{\|\cdot-\mathbf{t}\|}{h}\Big{)}$$

$$H_{0,+}^{\mathbf{t},h,\mathbf{v}}:{\int }_{{\mathbb{R}}^d}{\varphi }_{\mathbf{t},h}(\mathbf{b}){\partial }_{\mathbf{v}}{f}_{\mathbf{\beta }}(\mathbf{b})d\mathbf{b}\le 0$$

$$H_{0,-}^{\mathbf{t},h,\mathbf{v}}:{\int }_{{\mathbb{R}}^d}{\varphi }_{\mathbf{t},h}(\mathbf{b}){\partial }_{\mathbf{v}}{f}_{\mathbf{\beta }}(\mathbf{b})d\mathbf{b}\ge 0$$

$$H_{0,-}^{{\mathbf{t}}_j,h,{\mathbf{v}}_j}$$

$$f_{\mathbf{\beta }}(\mathbf{b})\ge c_{\mathbf{\beta }}\text{ for all }\mathbf{b}\in [{\mathbf{a}}_1,{\mathbf{a}}_2].$$

$$\tilde{\varphi }(z):={\int }_0^{\infty }{r}^{d-2}\frac{\partial }{\partial z}\varphi \left(\sqrt{z^2+r^2}\right)dr\text{mbox}withz\in \mathbb{R}.$$

$$\displaystyle\Lambda(\partial_{\mathbf{v}}\phi_{\mathbf{t},h})(s,{{\bm{\theta}}})=\frac{\langle{{\bm{\theta}}},\mathbf{v}\rangle}{h^{d+1}}(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})\Big{(}\frac{s-\langle\mathbf{t},{{\bm{\theta}}}\rangle}{h}\Big{)},$$

$$\displaystyle T_{\mathbf{t},h,\mathbf{v}}:=\frac{1}{n\sqrt{h}}\sum_{i=1}^{n}\frac{\langle{\bm{\Theta}}_{i},\mathbf{v}\rangle}{f_{{\bm{\Theta}}}({\bm{\Theta}}_{i})}(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})\Big{(}\frac{S_{i}-\langle\mathbf{t},{\bm{\Theta}}_{i}\rangle}{h}\Big{)}.$$

$$\displaystyle\mathbb{E}\big{[}T_{\mathbf{t},h,\mathbf{v}}\big{]}=\int_{\mathbb{S}^{d-1}}\int_{\mathbb{R}}h^{d+1/2}\Lambda(\partial_{\mathbf{v}}\phi_{\mathbf{t},h})(s,{{\bm{\theta}}})f_{S|{\bm{\Theta}}}(s|{{\bm{\theta}}})dsd{{\bm{\theta}}},$$

$$\displaystyle\mathbb{E}\big{[}T_{\mathbf{t},h,\mathbf{v}}\big{]}=-c_{d}h^{d+1/2}\int_{\mathbb{R}^{d}}\phi_{\mathbf{t},h}(\mathbf{b})\partial_{\mathbf{v}}f_{{\bm{\beta}}}(\mathbf{b})d\mathbf{b}.$$

$$\displaystyle\widehat{T}_{\mathbf{t},h,\mathbf{v}}:=\frac{1}{n\sqrt{h}}\sum_{i=1}^{n}\frac{\langle{\bm{\Theta}}_{i},\mathbf{v}\rangle}{\widetilde{f}_{{\bm{\Theta}}}({\bm{\Theta}}_{i})}(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})\Big{(}\frac{S_{i}-\langle\mathbf{t},{\bm{\Theta}}_{i}\rangle}{h}\Big{)}.$$

$$\displaystyle f_{\bm{\Theta}}({\bm{\theta}})=\frac{1}{2|\theta_{1}|^{d}}f_{\mathbf{X}}\Big{(}\frac{\theta_{2}}{\theta_{1}},\ldots,\frac{\theta_{d}}{\theta_{1}}\Big{)},$$

$$\displaystyle f_{\mathbf{X}}(\bm{x})=\frac{\Gamma(d/2)}{\pi^{d/2}|\Sigma|^{1/2}\big{(}1+(\bm{x}-\bm{\mu})^{\top}\Sigma^{-1}(\bm{x}-\bm{\mu})\big{)}^{d/2}}\quad\,\textnormal{for}\quad\bm{x}\in\mathbb{R}^{d-1},$$

$$\displaystyle f_{\bm{\Theta}}({\bm{\theta}})=\frac{\Gamma(d/2)}{2\pi^{d/2}|\Sigma|^{1/2}\big{(}\theta_{1}^{2}+((\mathrm{sgn}(\theta_{1})\theta_{j}-|\theta_{1}|\mu_{j})_{j=2}^{d})^{\top}\Sigma^{-1}((\mathrm{sgn}(\theta_{1})\theta_{j}-|\theta_{1}|\mu_{j})_{j=2}^{d})\big{)}^{d/2}},$$

$$\nu :\{\begin{array}{ll}\mathcal{B}(\mathcal{Z})& \to {\mathbb{R}}_0^{+},\\ E& \mapsto \nu (E)={\int }_{{\mathbb{S}}^{d-1}}{\int }_{\mathbb{R}}\text{mathds}{1}_E(\mathbf{\theta },s)dsd\mathbf{\theta }.\end{array}$$

$$W(E_1)\sim \mathcal{N}(0,\nu (E_1)),W(E_1\cup E_2)=W(E_1)+W(E_2)\text{ a.s.}\text{and}W(E_1)\perp W(E_2)$$

$$\displaystyle\int_{E}\,W(dsd{\bm{\theta}})=W(E)\sim\mathcal{N}(0,\nu(E))\quad\text{for}\quad E\in\bigl{(}\mathcal{B}(\mathcal{Z})\bigr{)}_{\nu},$$

$$W(f):={\int }_{{\mathbb{S}}^{d-1}}{\int }_{\mathbb{R}}f(s,\mathbf{\theta })W(dsd\mathbf{\theta })\sim \mathcal{N}(0,\Vert f{\Vert }_{L^2(\nu )})\text{for}f\in L^2(\nu ),$$

$$\displaystyle X_{\mathbf{t},h,\mathbf{v}}=h^{-1/2}\int_{\mathbb{S}^{d-1}}\int_{\mathbb{R}}\langle{\bm{\theta}},\mathbf{v}\rangle(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})\Big{(}\frac{s-\langle\mathbf{t},{\bm{\theta}}\rangle}{h}\Big{)}\frac{\sqrt{f_{S,{\bm{\Theta}}}(s,{{\bm{\theta}}})}}{{f_{{\bm{\Theta}}}({{\bm{\theta}}})}}W(dsd{\bm{\theta}}).$$

Full text

Keywords:

Gaussian approximation; mode detection; monotonicity; multiscale statistics; shape constraints; Radon transform; ill-posed problems.

1 Introduction

In the random coefficients model, $n$ i.i.d. random vectors $(\mathbf{X}_{i},Y_{i}),$ $i=1,\ldots,n$ are observed, with $\mathbf{X}_{i}=(X_{i,1},\ldots,X_{i,d})$ a $d$-dimensional vector of design variables and

[TABLE]

The unobserved random coefficients ${\bm{\beta}}_{i}=(\beta_{i,1},\ldots,\beta_{i,d}),$ $i=1,\ldots,n,$ are i.i.d. realizations of an unknown $d$-dimensional distribution $F_{\bm{\beta}}$ with Lebesgue density $f_{\bm{\beta}}.$ Design variables and random coefficients are assumed to be independent. The statistical task is to recover properties of the joint density $f_{\bm{\beta}},$ which is assumed to belong to some nonparametric class. In this work, we derive tests for increases and modes of $f_{\bm{\beta}}.$

For $d=1,$ the random coefficients model simplifies to nonparametric density estimation. For $d>1,$ recovery of $f_{\bm{\beta}}$ is an inverse problem with ill-posedness depending on the distribution of the design vectors $\mathbf{X}_{i}.$ If the design is sufficiently regular, the inverse problem is mildly ill-posed. Otherwise, the model can be severely ill-posed or even be non-identifiable. In this work, we study the mildly ill-posed regime and consider in particular the random coefficients model with random intercept

[TABLE]

which can be obtained from (1.1) setting $X_{i,1}=1,$ almost surely.

Random coefficients models appear in econometrics and epidemiology and are used to model unobserved heterogeneity in the population. While the standard linear regression model accounts for unobserved heterogeneity only by an intercept that varies across the population, the random coefficients model allows in addition that different individuals have different slopes. Applications in epidemiology are considered by Greenland, (2000); Gustafson and Greenland, (2006). In economics, random coefficients models are frequently used to evaluate panel data, cf. Hsiao, (2014) or Hsiao and Pesaran, (2004), Chapter 6, for an overview. Modeling and estimating consumer demand in industrial organization and marketing often makes use of random coefficients Berry et al., (1995); Petrin, (2002); Nevo, (2001); Berry and Pakes, (2007); Dubé et al., (2012). In all these works, parametric assumptions on $f_{{\bm{\beta}}}$ are imposed. Recently, nonparametric approaches for random coefficients became popular in microeconometrics Hoderlein et al., (2010); Masten, (2017); Hoderlein et al., (2015); Dunker et al., (2017), frequently combined with binary choice Ichimura and Thompson, (1998); Gautier and Hoderlein, (2012); Gautier and Kitamura, (2013); Masten and Torgovitsky, (2014); Dunker et al., (2013); Fox and Gandhi, (2016); Dunker et al., (2018), among others.

The random coefficients model also includes quantum homodyne tomography. In this case, we observe an angle $\Phi_{i}$ and

[TABLE]

with $(Q_{i},P_{i})$ i.i.d. random variables which are unobserved and independent of $\Phi_{i}.$ The angles $\Phi_{i}$ can be chosen by the experimenter and are typically uniform on $[0,\pi].$ The interest is in reconstruction of the Wigner function which takes the role of the joint density of $(Q_{i},P_{i}).$ Because $P_{i}$ and $Q_{i}$ are not jointly observable, the Wigner function can take negative values. For more on quantum homodyne tomography and the Wigner function, see Butucea et al., (2007).

We propose a nonparametric test for shape information of the joint density $f_{\bm{\beta}}$ in the random coefficients model. The focus will be on a test for directional derivatives and modes. The nonparametric estimation theory for $f_{{\bm{\beta}}}$ has been developed in Beran and Hall, (1992); Beran et al., (1996); Feuerverger and Vardi, (2000); Hoderlein et al., (2010). Due to the ill-posedness of the problem and the curse of dimensionality induced by $d,$ pointwise estimation rates are slow. The reason is that small perturbations in the signal are indistinguishable given the data. Nevertheless, we can get good detection rates for larger features, such as an accentuated mode or a strong increase in the joint density along some direction. From a practical point of view, the relevant information regarding an unknown density is typically its shape rater than its precise, full reconstruction. It is therefore essential to recover increases/decreases and the modes of a density. If, say, two modes in the joint density of two random quantities are detected, this indicates that two different groups can be identified. Hence, shape information allows to interpret a given dataset.

Larger features of the density will also be discovered by a nonparametric estimator even if it suffers from slow pointwise convergence. There are, however, two important reasons why a testing approach might be more appropriate. Firstly, with a significance test of level $\alpha$ we can conclude that with probability $1-\alpha$ a detected feature is not an artifact. Secondly, for an estimator we need to pick one bandwidth or smoothing parameter while detection of different features might require different bandwidth choices depending on the size of the hidden features themselves. Indeed, a short and steep increase will be best detected on a small scale whereas for finding a longer and less strong increase the choice of a larger bandwidth is beneficial. Using multiple testing methods, it is possible to combine a whole range of smoothness parameters into one test and to adapt to different shapes of features.

We construct a so called multiscale test, aggregating single tests on different scales and directions. Multiscale tests can be viewed as a multiple testing procedure specifically designed for nonparametric models. Given a model, the theoretical challenge is to prove that a multiscale statistic can be approximated by a distribution free statistic which is independent of the observations. This allows us then to compute quantiles and to find approximations for the critical values of the multiscale statistic. So far, qualitative feature detection based on multiscale statistics has been studied for various nonparametric models, including the Gaussian white noise model Dümbgen and Spokoiny, (2001), density estimation Dümbgen and Walther, (2008) and deconvolution Schmidt-Hieber et al., (2013). In multivariate settings the classical KMT approximation suffers from the curse of dimensionality which then leads to very restrictive conditions on the usable scales. Instead, very recent results on Gaussian approximations of suprema of multivariate empirical processes developed by Chernozhukov et al., (2017) can be used (Eckle et al., 2017a, ; Proksch et al.,, 2016, in the context of multivariate deconvolution and multivariate linear inverse problems with additive noise, respectively). In this work, we extend these techniques. The main difficulties are twofold. First, we need to derive specific properties of the inverse Radon transform for general dimension $d$. Second, in contrast to the other works on multiscale inference, no distribution free approximation can be obtained and we therefore need to study the approximating process if several unobserved functions are replaced by estimators.

In order to study the power of the multiscale test, a theoretical detection bound and numerical simulations are provided. The theoretical result gives conditions under which a mode can be detected. In a numerical simulation study, we investigate the power of the test for increases/decreases along some direction and mode detection in dependence on the sample size and the design variables. We also analyze real consumer demand data from the British Family Expenditure Survey.

Let us briefly summarize related literature on testing in the random coefficients model. Under a parametric assumption on the density $f_{\bm{\beta}},$ Beran, (1993) considers goodness-of-fit testing and Swamy, (1970); Andrews, (2001) test whether some of the random coefficients are deterministic. The only test based on a nonparametric assumption was proposed recently by Breunig and Hoderlein, (2018). It allows to assess whether a given set of data follows the random coefficients model.

This paper is organized as follows. In Section 2, we describe the connection between the random coefficients model and the Radon transform. Rewriting the model as an inverse problem in terms of the Radon transform reveals the ill-posed nature of the model. This allows us to construct and to analyze the multiscale test in Section 3. In this part we also derive the asymptotic theory of the estimator and obtain theoretical detection bounds. In Section 4 the test is studied for simulated data. As a real data example, consumer demand is analyzed in Section 5. Proofs and technicalities are deferred to a supplement. An R package and Python code is available as a supplement as well, https://arxiv.org/abs/1704.01066.

Notation: Throughout the paper, vectors are displayed by bold letters, e.g. $\mathbf{X},{\bm{\beta}}.$ Inequalities between vectors are understood componentwise. The Euclidean norm on $\mathbb{R}^{d}$ is denoted by $\|\cdot\|$ and the corresponding standard inner product by $\langle\cdot,\cdot\rangle.$ We further denote by $\mathbf{e}_{1},\ldots,\mathbf{e}_{d}\in\mathbb{R}^{d}$ the standard ON-basis of the $d$-dimensional Euclidean space, $\mathbb{S}^{d-1}$ denotes the unit sphere in $\mathbb{R}^{d}$ and we write $\mathcal{Z}$ for the cylinder $\mathcal{Z}=\mathbb{R}\times\mathbb{S}^{d-1}$. Furthermore, we write $\mathbf{v}$ for any direction $\mathbf{v}=\sum_{j=1}^{d}\mathrm{v}_{j}\mathbf{e}_{j}\in\mathbb{S}^{d-1}$. For two positive sequences $(a_{n})_{n},$ $(b_{n})_{n},$ $a_{n}\lesssim b_{n}$ or $b_{n}\gtrsim a_{n}$ mean that for some positive constant $C,$ $a_{n}\leq Cb_{n}$ for all $n.$ As usual, we write $a_{n}\asymp b_{n}$ if $a_{n}\lesssim b_{n}$ and $b_{n}\lesssim a_{n}.$

2 The random coefficients model as an inverse problem

The random coefficients model can be written in terms of the Radon transform (cf. Beran et al.,, 1996). This allows us then to interpret the model as an inverse problem. In this section, we summarize the main steps and review relevant results on the inversion of the Radon transform. Let $H^{s}$ denote the $L^{2}$-Sobolev space. The Radon transform is the operator $R:H^{s}(\mathbb{R}^{d})\rightarrow H^{s+\frac{d-1}{2}}(\mathcal{Z}),$ with

[TABLE]

and $\mu_{d-1}$ the surface measure on the $(d-1)$-dimensional hyperplane $\{\mathbf{b}\in\mathbb{R}^{d}:\langle\mathbf{b},{\bm{\theta}}\rangle=s\}$. The Radon transform maps therefore a function to all its integrals over hyperplanes parametrized by $(s,{\bm{\theta}})\in\mathcal{Z}.$ The figure above shows the parametrization in two dimensions.

For the connection between the Radon transform and the random coefficients model (1.1) consider the normalized observations

[TABLE]

The random vectors ${\bm{\Theta}}_{i}$ take values in the $(d-1)$-dimensional sphere $\mathbb{S}^{d-1}$. In the random coefficients model with intercept (1.2), ${\bm{\Theta}}_{i}$ is always in the upper hemisphere, i.e. the first component of ${\bm{\Theta}}_{i}$ is positive. In this case, we extend the distribution of ${\bm{\Theta}}_{i}$ to the whole sphere by randomizing the signs of the design variables. For this purpose, we generate independent random variables $\zeta_{i}$, $i=1,\ldots,n$, with $\mathbb{P}(\zeta_{i}=1)=\mathbb{P}(\zeta_{i}=-1)=1/2,$ which are independent of the data $(\mathbf{X}_{i},Y_{i}),$ $i=1,\ldots,n$, and define $S_{i}:=\zeta_{i}Y_{i}/\|\mathbf{X}_{i}\|$ and ${\bm{\Theta}}_{i}:=\zeta_{i}\mathbf{X}_{i}/\|\mathbf{X}_{i}\|.$ Independent of the symmetrization, we have $S_{i}=\langle{\bm{\Theta}}_{i},{\bm{\beta}}_{i}\rangle$. The conditional distribution of $S_{1}|{\bm{\Theta}}_{1}$ is therefore

[TABLE]

and the conditional density becomes

[TABLE]

Recall that we have access to an i.i.d. sample $(S_{i},{\bm{\Theta}}_{i})$ of the joint density $f_{S,{\bm{\Theta}}}$. This allows for nonparametric estimation of $f_{S|{\bm{\Theta}}}$ via $f_{S|{\bm{\Theta}}}(s|{\bm{\theta}})=f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})/f_{\bm{\Theta}}({\bm{\theta}})$. Applying the inverse Radon transform to this estimate gives an estimator for the joint density $f_{\bm{\beta}}.$ This inversion scheme suffers from two sources of ill-posedness. Firstly, dividing by $f_{\bm{\Theta}}$ might result in very unstable reconstructions if $f_{\bm{\Theta}}$ is small. This happens if the normalized design variables ${\bm{\Theta}}_{i}=\mathbf{X}_{i}/\|\mathbf{X}_{i}\|$ systematically miss observations from some directions. In this case the problem becomes unevenly harder and only logarithmic convergence rates can be obtained (see Davison,, 1983; Frikel,, 2013; Hohmann and Holzmann,, 2016). When the support of ${\bm{\Theta}}_{i}$ does not contain an open ball, $f_{\bm{\beta}}$ might be non-identifiable. Secondly, even with regularity on the distribution of the design, the Radon inversion is known to be an ill-posed problem with degree of ill-posedness $(d-1)/2$. Hence, regularization of the inversion scheme is necessary.

In this work, we study the mildly ill-posed case where the random directions ${\bm{\Theta}}_{i},$ $i=1,\ldots,n,$ are sufficiently regularly distributed over the sphere and the ill-posedness is only due to the inversion of the Radon transform. The precise assumptions on the design are stated in Section 3.2.

Our approach makes use of the following explicit inversion formula of the Radon transform. Define the operator $\Lambda$ via

[TABLE]

where $\mathcal{H}_{d}$ denotes the identity for $d$ odd and the Hilbert transform

[TABLE]

for $d$ even. Let $c_{d}^{-1}=(-1)^{({d-1})/{2}}2^{-d}\pi^{{1-d}}$ for $d$ odd and ${c}_{d}^{-1}=-(-1)^{{d}/{2}}2^{-d}\pi^{{1-d}}$ for $d$ even. If $\varphi$ is a Schwartz function on $\mathbb{R}^{d}$, then we have the inversion formula

[TABLE]

cf. Theorem 3.8 in Helgason, (2011). The so called back projection operator $R^{*}$ is the adjoint of the Radon transform with respect to the $L^{2}$ scalar product. Notice that our constant $c_{d}$ differs from the constant in Helgason, (2011) as we use the standard definition of the Hilbert transform and define $R^{*}$ as the adjoint of the Radon transform (as opposed to the dual transform).

3 Multiscale tests for qualitative features

3.1 Multiscale inference

The goal of this work is to derive confidence statements for qualitative features of the joint density of the random coefficients. In particular, we are interested in the detection of modes (local maxima) of the density. Following the approach of Schmidt-Hieber et al., (2013), we express the features in terms of differential operators. To be precise, for a collection of compactly supported, non-negative and sufficiently smooth test functions $\phi_{\mathbf{t},h}$ consider the integral

[TABLE]

for some directional vector $\mathbf{v}=(\rm v_{1},\ldots,\rm v_{d})^{\top}$ and $d\geq 2$. Since there should not be any favored direction, we consider in the following radially symmetric test functions,

[TABLE]

with a non-negative and sufficiently smooth kernel $\phi:[0,\infty)\rightarrow[0,\infty)$ with $\int_{0}^{\infty}\phi(u)du=1$ and support on $[0,1]$. Moreover, $\operatorname{Vol}(\mathbb{S}^{d-2})$ denotes the volume of the sphere $\mathbb{S}^{d-2}\subset\mathbb{R}^{d-1}$ and $\operatorname{Vol}(\mathbb{S}^{0}):=2$. Notice that $\phi_{\mathbf{t},h}$ is supported on the ball $B_{h}(\mathbf{t})$ with center $\mathbf{t}$ and radius $h.$ The normalization for $\phi_{\mathbf{t},h}$ turns out to be convenient but does not entail that $\phi_{\mathbf{t},h}$ integrates to one.

If the integral (3.1) is positive, there exists a subset of $B_{h}(\mathbf{t})$ with positive Lebesgue measure on which $\partial_{\mathbf{v}}f_{{\bm{\beta}}}$ is positive. On this subset, $f_{{\bm{\beta}}}$ is thus strictly increasing in direction $\mathbf{v}.$ Similarly, we can recover a decrease if the integral (3.1) is negative. To construct a statistical test for increases and decreases it is therefore natural to use an empirical counterpart of the functional defined in (3.1).

Let $\mathcal{T}=\{(\mathbf{t},h,\mathbf{v}):h\in(0,1],\mathbf{a}_{1}+h\leq\mathbf{t}\leq\mathbf{a}_{2}-h,\mathbf{v}\in\mathbb{S}^{d-1}\}$ for $\mathbf{a}_{1}\leq\mathbf{a}_{2}\in\mathbb{R}^{d}$, where the inequalities for the vectors $\mathbf{a}_{1},\mathbf{a}_{2},\mathbf{t}$ are understood componentwise. For statistical inference regarding the sign of the directional derivatives of $f_{\bm{\beta}}$, we fix a subset $\mathcal{T}_{n}\subset\mathcal{T}$ and test for all $(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}$ simultaneously the corresponding hypotheses of the form

[TABLE]

and

[TABLE]

For constructing global tests, we can now argue as in Eckle et al., 2017a . Our main interest are the following three global testing problems (i)-(iii).

(i) Testing for the presence of a mode at a fixed location. Tests for the hypotheses (3.3) and (3.4) can be used for the detection of specific shape constraints such as a mode at a given point $\mathbf{b}_{0}\in\mathbb{R}^{d}$. For this purpose, we consider several bandwidths/scales $h$ and for each $h$ consider pairs $(\mathbf{t}_{1},\mathbf{v}_{1}),\ldots,(\mathbf{t}_{p},\mathbf{v}_{p})$, where $\mathbf{v}_{j},$ $j=1,\ldots,p$, are directional vectors and the test locations $\mathbf{t}_{j}$ are points on the line $\{\mathbf{b}_{0}+r\mathbf{v}_{j}:r\geq h\}$ $(j=1,\ldots,p)$ in a neighborhood of $\mathbf{b}_{0}$. Inference for the presence of a mode at the point $\mathbf{b}_{0}$ can now be conducted by studying the testing problem

[TABLE]

with $h$ ranging over all chosen scales and $j=1,\ldots,p.$ Level and power of the mode test (3.5) for different designs are reported in Section 4.2. In Section 4.3, we also show that it is essential to include several bandwidths/scales $h$ in order to separate modes which are close.

(ii) A global testing procedure for all modes. Simultaneous tests for the hypotheses (3.3) and (3.4) can be used for a global testing procedure to detect all modes of the density on a domain. Compared to the previous case, we search for evidence over a range of different $\mathbf{b}_{0}$ which then inflates the number of local tests.

(iii) A graphical representation of the local monotonicity behavior for bivariate densities. Let $d=2$ and define a subset $\mathcal{T}_{n}=\{(\widetilde{\mathbf{t}}_{j},h_{0},\widetilde{\mathbf{v}}_{j}):j=1,\ldots,p\}$ for a fixed scale $h_{0}$ of the form $\mathcal{T}_{n}=\mathcal{T}_{\mathbf{t}}\times\{h_{0}\}\times\mathcal{T}_{\mathbf{v}}$, where $\mathcal{T}_{\mathbf{t}}$ contains the $p/|\mathcal{T}_{\mathbf{v}}|$ vertices of an equidistant grid of width $2h_{0}$ and $\mathcal{T}_{\mathbf{v}}$ contains the directions. We restrict the testing procedure to one fixed scale for an easy to read graphical representation. We consider four equidistant directions on $\mathbb{S}^{1}$ given by $\mathcal{T}_{\mathbf{v}}=\{\mathbf{v}_{1},-\mathbf{v}_{1},\mathbf{v}_{2},-\mathbf{v}_{2}\}.$ Since $\mathcal{T}_{\mathbf{v}}=-\mathcal{T}_{\mathbf{v}}$, we have symmetry in the hypotheses, i.e. $H_{0,+}^{\widetilde{\mathbf{t}}_{j},h_{0},\widetilde{\mathbf{v}}_{j}}=H_{0,-}^{\widetilde{\mathbf{t}}_{j},h_{0},-\widetilde{\mathbf{v}}_{j}}$. Therefore, we test only $H_{0,-}^{\widetilde{\mathbf{t}}_{j},h_{0},\widetilde{\mathbf{v}}_{j}}$ for all triples $(\widetilde{\mathbf{t}}_{j},h_{0},\widetilde{\mathbf{v}}_{j})\in\mathcal{T}_{n}.$ Figure 1 displays an example for the test outcome with the hypotheses in (3.3) and $\mathcal{T}_{n}$ as above. An arrow in a direction $\widetilde{\mathbf{v}}_{j}$ at a location $\widetilde{\mathbf{t}}_{j}$ represents a rejection of the corresponding hypothesis $H_{0,-}^{\widetilde{\mathbf{t}}_{j},h_{0},\widetilde{\mathbf{v}}_{j}}$ and provides an indication of a negative directional derivative of $f_{\bm{\beta}}$ in direction $\widetilde{\mathbf{v}}_{j}$ at the location $\widetilde{\mathbf{t}}_{j}$. Thus, Figure 1 provides strong evidence that the density is trimodal with modes close to the locations $(-0.5,-0.5)^{\top}$, $(1.5,-0.5)^{\top}$, and $(0.5,1.5)^{\top}$. A detailed description of the settings used to generate Figure 1 and an analysis of the results is given in Section 4.1.

We now derive an empirical counterpart of the functional (3.1) in terms of the Radon transform $Rf_{\bm{\beta}}$. We make the following assumptions for the inversion of the Radon transform.

Assumption 1**.**

Suppose that the function $\phi$ in (3.2) is $(d+2)$-times continuously differentiable with $\phi^{\prime}(0)=\phi^{\prime\prime}(0)=0$.

Assumption 2**.**

Suppose that the density $f_{\bm{\beta}}$ is compactly supported, continuously differentiable and bounded from below in the test region by a constant $c_{\bm{\beta}}>0$

[TABLE]

Assumption 2 is too restrictive for quantum homodyne tomography (model (1.3)), where the density $f_{\bm{\beta}}$ is given by the Wigner function. The Wigner function can take negative values and is not compactly supported. In this case, we replace Assumption 2 by the following conditions.

Assumption 2′.**

Suppose that $f_{\bm{\beta}}$ is continuously differentiable and, for some $\gamma^{\prime},\varepsilon>0$,

1. (i)

$\mathbf{b}\mapsto\|\mathbf{b}\|^{d+\varepsilon}\displaystyle|f_{\bm{\beta}}(\mathbf{b})|$* is bounded;* 2. (ii)

$\displaystyle|f_{\bm{\beta}}(\mathbf{b})-f_{\bm{\beta}}(\mathbf{b}^{\prime})|\lesssim\frac{\|\mathbf{b}-\mathbf{b}^{\prime}\|^{\gamma^{\prime}}}{(1+\min\{\|\mathbf{b}\|,\|\mathbf{b}^{\prime}\|\})^{d+\gamma^{\prime}+\varepsilon}}$* for all $\mathbf{b},\mathbf{b}^{\prime}\in\mathbb{R}^{d}$;* 3. (iii)

There exist constants $\delta,c_{\bm{\beta}}>0$ such that for every hyperplane $P\subset\mathbb{R}^{d}$ with $P\cap[\mathbf{a}_{1}-\delta,\mathbf{a}_{2}+\delta]\neq\emptyset$, $\int_{P}f_{\bm{\beta}}(\mathbf{b})d\mu_{d-1}(\mathbf{b})\geq c_{\bm{\beta}}.$

Under the assumptions above the inversion formula (2.3) holds for partial derivatives of the test functions $\partial_{\mathbf{v}}\phi_{\mathbf{t},h}$. This is a direct consequence of Theorem 3.8 in Helgason, (2011). The following lemma analyzes the structure of a partial derivative of the test function transformed by the operator $\Lambda$ introduced in (2.2) and how this transform depends on $h$.

Lemma 3.1**.**

Work under Assumption 1 and let

[TABLE]

Then

[TABLE]

where $\Lambda$ and $\phi_{\mathbf{t},h}$ are defined in (2.2) and (3.1), respectively. Moreover,

[TABLE]

For a given triple $(\mathbf{t},h,\mathbf{v})\in\mathcal{T}$ we study the statistic

[TABLE]

By Lemma 3.1, the expectation of this statistic can be written as

[TABLE]

with $d{\bm{\theta}}$ being the surface measure on $\mathbb{S}^{d-1}$, i.e. $|S_{\bm{\Theta}}|=\int_{S_{\bm{\Theta}}}\,d{\bm{\theta}}$ for any measurable $S_{\bm{\Theta}}\subseteq\mathbb{S}^{d-1}.$ By an application of the inversion formula introduced in (2.3) and Lemma 5.1 in Helgason, (2011), we obtain

[TABLE]

Up to rescaling, $T_{\mathbf{t},h,\mathbf{v}}$ is thus an empirical counterpart of the functional defined in (3.1).

The statistic $T_{\mathbf{t},h,\mathbf{v}}$ depends on the density $f_{\bm{\Theta}}.$ In quantum homodyne tomography this density is known. For many other applications, however, $f_{\bm{\Theta}}$ needs to be estimated from the data. In this case, we use a standard cut-off kernel density estimator $\widetilde{f}_{{\bm{\Theta}}}$ for $f_{\bm{\Theta}}$ based on an additional sample $(S_{i},{\bm{\Theta}}_{i}),$ $i=n+1,\ldots,2n$ which is independent of $(S_{i},{\bm{\Theta}}_{i}),$ $i=1,\ldots,n$. See also Appendix A for the definition of the estimator. We replace $f_{\bm{\Theta}}$ by its estimator and consider the test statistic

[TABLE]

3.2 Assumptions on the design

As mentioned in Section 2, the inverse problem might become severely ill-posed or non-identifiable if the density $f_{\bm{\Theta}}$ approaches zero for some directions. This section provides conditions on the design which ensure that $f_{\bm{\Theta}}$ has positive Hölder smoothness and is bounded from below and above. These results are of independent interest.

In the random coefficients model (1.1), the density $f_{\bm{\Theta}}$ can be expressed in terms of the density $f_{\mathbf{X}}$ via $f_{\bm{\Theta}}({\bm{\theta}})=\int_{0}^{\infty}r^{d-1}f_{\mathbf{X}}(r{\bm{\theta}})dr$. To enforce that $f_{\bm{\Theta}}$ is bounded from below we restrict ourselves to designs where $\int_{0}^{\infty}r^{d-1}f_{\mathbf{X}}(r{\bm{\theta}})dr$ is bounded away from [math]. The formula also allows to relate the smoothness of $f_{\mathbf{X}}$ to the smoothness of $f_{\bm{\Theta}}.$

Although, the random coefficients model with intercept (1.2) could be viewed as a special case of the more general model (1.1) it requires a different set of assumptions. For model (1.2), we write $f_{\mathbf{X}}$ as a function of $\bm{x}=(x_{2},\ldots,x_{d})\in\mathbb{R}^{d-1}$ and obtain

[TABLE]

see Appendix B for a proof. A necessary condition to ensure that $\inf_{{\bm{\theta}}}f_{\bm{\Theta}}({\bm{\theta}})>0,$ is given by $f_{\mathbf{X}}(\bm{x})\gtrsim\|\bm{x}\|^{-d}$ as $\|\bm{x}\|\rightarrow\infty$. This corresponds to Cauchy-type tails of the design variables. Thinner tails will increase the ill-posedness of the problem. In order to avoid very technical proofs, we consider in the random coefficients model with intercept only the case where $(X_{i,2},\ldots,X_{i,d})$ follows a multivariate Cauchy distribution, i. e.,

[TABLE]

with $\bm{\mu}\in\mathbb{R}^{d-1}$ and $\Sigma\in\mathbb{R}^{(d-1)\times(d-1)}$ a symmetric and positive definite matrix. We can compute $f_{\bm{\Theta}}$ explicitly using (3.9)

[TABLE]

where $\mathrm{sgn}(\cdot)$ denotes the signum function. In this case, $f_{\bm{\Theta}}$ is bounded from above and below and is continuously differentiable on the hemispheres $\mathbb{S}^{d-1}_{+}:=\{{\bm{\theta}}\in\mathbb{S}^{d-1}\mid\theta_{1}>0\}$ and $\mathbb{S}^{d-1}_{-}:=\{{\bm{\theta}}\in\mathbb{S}^{d-1}\mid\theta_{1}<0\}$. In particular, if $(X_{i,2},\ldots,X_{i,d})$ is standard Cauchy, then the density $f_{\bm{\Theta}}$ is constant. This leads to the following assumptions on the design.

Assumption 3**.**

In model (1.1), suppose that

[TABLE]

In model (1.2), assume that $f_{\mathbf{X}}$ is of the form (3.10) with $\Sigma$ a symmetric and positive definite matrix.

In quantum homodyne tomography we set a global $\gamma$ equal to the minimum of $\gamma^{\prime}$ from Assumption 2’ and $\gamma$ from Assumption 3.

It is important to notice that statistical testing in the random coefficients model relies on two unrelated sets of assumptions. Firstly, there are assumptions on the density of the random coefficients as introduced in Section 3.1. Restrictions of this kind are common in statistical inference for an unknown density. On the other hand, there are assumptions (see Assumption 3 above) on the design. These assumptions control the ill-posedness of the problem.

Note that Assumption 3 (ii) can be weakened to $\int_{-\infty}^{\infty}r^{d-1}f_{\mathbf{X}}(r{\bm{\theta}})dr\geq c>0\quad\text{for all }{\bm{\theta}}\in\mathbb{S}^{d-1}$. For example, if the support of $\Theta$ is a hemisphere this condition can hold while Assumption 3 (ii) is violated. This relaxation can be achieved by multiplying independently generated $\zeta_{i}$ to $\Theta_{i}$ and $S_{i}$ as proposed for model (1.2) in Section 2.

3.3 Asymptotic properties

This section presents the main theoretical result of the paper stating that the standardized and properly calibrated test statistic (3.8) can be uniformly approximated by a maximum of a Gaussian process. For that we need the definition of a Gaussian process on the cylinder $\mathcal{Z}$. To this end, let $\mathcal{B}(\mathcal{Z})$ be the Borel $\sigma$-algebra on $\mathcal{Z}$. Define the $\sigma$-finite measure

[TABLE]

Let \bigl{(}\mathcal{B}(\mathcal{Z})\bigr{)}_{\nu} denote the collection of all sets of finite $\nu$-measure and let $W$ denote Gaussian $\nu$-noise on \bigl{(}\mathcal{B}(\mathcal{Z}),\nu\bigr{)}. For disjoint sets E_{1},E_{2}\in\bigl{(}\mathcal{B}(\mathcal{Z})\bigr{)}_{\nu} this implies

[TABLE]

(Adler and Taylor,, 2007, Chapter 1.4.3). $W$ is a random, finitely additive, signed measure. Integration w.r.t. $W$ can be defined similarly to Lebesgue-integration, starting with a definition for simple functions and an extension to general $f\in L^{2}(\nu)$ via approximation by simple functions in the $L^{2}$-limit. Integration with respect to $W$ yields

[TABLE]

[TABLE]

and $\mathrm{Cov}(W(f),W(g))=\langle W(f),W(g)\rangle_{L^{2}(\mathbb{P})}=\langle f,g\rangle_{L^{2}(\nu)}$ for $f,g\in L^{2}(\nu),$ where $L^{k}(\mathbb{P})$ denotes the collection of all random variables whose first $k$ absolute moments exist. For more details, cf. Adler and Taylor, (2007), Chapter 5.2.

Let us provide some heuristic for the Gaussian approximation of $T_{\mathbf{t},h,\mathbf{v}}$. The process $(\mathbf{t},h,\mathbf{v})\mapsto\sqrt{n}(T_{\mathbf{t},h,\mathbf{v}}-\mathbb{E}[T_{\mathbf{t},h,\mathbf{v}}])$ has in the important case $\mathbb{E}[T_{\mathbf{t},h,\mathbf{v}}]=0$ the same mean and covariance structure as the Gaussian process

[TABLE]

In the proof of Theorem 3.2 below we show that the expectation $\mathbb{E}[T_{\mathbf{t},h,\mathbf{v}}]$ is asymptotically negligible in the limit process. The test statistic and the Gaussian process depend, however, on the unknown densities $f_{S,{\bm{\Theta}}}$ and $f_{{\bm{\Theta}}}$ which have to be estimated from the second part of the sample. To this end, we use the standard cut-off kernel density estimates $\widetilde{f}_{{\bm{\Theta}}}$ and $\widetilde{f}_{S,{\bm{\Theta}}}$ defined in Appendix A. For Gaussian $\nu$-noise $W$ that is independent of the data let

[TABLE]

and

[TABLE]

The Gaussian approximation result for the family of test statistics $\widehat{T}_{\mathbf{t},h,\mathbf{v}}$ holds for a finite subset $\mathcal{T}_{n}\subset\mathcal{T}$. Its cardinality may, however, grow polynomially of arbitrary degree with the sample size. Moreover, the range of bandwidths must be bounded from above and below by $h_{\max}$ and $h_{\min}$, both converging to zero as $n$ goes to infinity. The precise conditions are summarized in the following assumption.

Assumption 4**.**

When working in model (1.1), let $\gamma$ be as defined in Assumption 3. When considering model (1.2), let $\gamma=1.$ Set $h_{\min}:=\min\{h:(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}\}$ and $h_{\max}:=\max\{h:(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}\}.$ Suppose that $|\mathcal{T}_{n}|=p\lesssim n^{L}$ for some $L>0$ and $h_{\,\textnormal{max}}\lesssim\log(n)^{-14\gamma/(d-1)-5}n^{2\gamma/(d-1)-1}\wedge o(1)$, $h_{\min}\gtrsim n^{-1+\varepsilon}$ for some $\varepsilon>0$.

Let $\mathcal{A}_{p}$ be the set of half-open hyperrectangles in $\mathbb{R}^{p}$, i.e. every $A\in\mathcal{A}_{p}$ has the representation $A=\{\bm{x}\in\mathbb{R}^{p}:-\infty<\bm{x}\leq\bm{a}\}$ for some $\bm{a}\in\mathbb{R}^{p}$. For finite sets $S_{n}$ and two stochastic processes $(X_{s,n})_{s\in S_{n}}$ and $(\widetilde{X}_{s,n})_{s\in S_{n}},$ which are defined on the same probability space, we write

[TABLE]

if $\lim_{n}\sup_{A\in\mathcal{A}_{|S_{n}|}}|\mathbb{P}((X_{s,n})_{s\in S_{n}}\in A)-\mathbb{P}((\widetilde{X}_{s,n})_{s\in S_{n}}\in A)|=0.$

Theorem 3.2**.**

For the calibration of the standardized statistic define

[TABLE]

Then, under Assumptions 1-4,

[TABLE]

Furthermore, almost surely conditionally on the estimators $\widetilde{f}_{S,{\bm{\Theta}}}$ and $\widetilde{f}_{\bm{\Theta}},$ the limit distribution is bounded in probability by a constant that is independent of the sample size $n$.

3.4 Construction of the multiscale test

With the previous theorem, we can now construct simultaneous statistical tests for the hypotheses (3.3) and (3.4). If the constant $c_{d}$ is positive then the method consists of rejecting the hypotheses $H_{0,+}^{\mathbf{t},h,\mathbf{v}}$ in (3.3) for small values of $\widehat{T}_{\mathbf{t},h,\mathbf{v}}$ and rejecting $H_{0,-}^{\mathbf{t},h,\mathbf{v}}$ in (3.4) for large values of $\widehat{T}_{\mathbf{t},h,\mathbf{v}}$, and vice versa if $c_{d}$ is negative. Theorem 3.2 is used to control the multiple level of the tests. Let $\alpha\in(0,1)$ and denote by $\kappa_{n}(\alpha)$ the smallest number such that

[TABLE]

By Theorem 3.2, $\kappa_{n}(\alpha)$ is bounded uniformly with respect to $n$. Define for $(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}$ the quantiles

[TABLE]

and reject the hypothesis (3.3), if

[TABLE]

Similarly, the hypothesis (3.4) is rejected, whenever

[TABLE]

Theorem 3.3**.**

Let the assumptions of Theorem 3.2 hold and assume that the tests (3.15) and (3.16) for the hypotheses (3.3) and (3.4) are performed simultaneously for all $(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}$. The probability of at least one false rejection of any of the tests is asymptotically at most $\alpha$, i.e.

[TABLE]

for $n\rightarrow\infty$.

Based on the previous result, we now propose a method for the detection and localization of modes on a subdomain. For convenience, we only consider the case of a hyperrectangle $[\mathbf{a}_{1},\mathbf{a}_{2}]\subset\mathbb{R}^{d}$. We study the case that there is a mode $\mathbf{b}_{0}$ in the interior $(\mathbf{a}_{1},\mathbf{a}_{2}).$ For the multiscale test to have power we need that the set of local tests is rich enough. This can be expressed in terms of conditions on $\mathcal{T}_{n}.$ Let $H(\mathcal{T}_{n})$ be the set of all scales/bandwidths such that for every scale $h$ in this set all triples $(\mathbf{t},h,\mathbf{v})$ are in $\mathcal{T}_{n},$ where $\mathbf{t}$ ranges over all grid points of an equidistant grid with component wise mesh size $h$ in the hyperrectangle $[\mathbf{a}_{1}+h,\mathbf{a}_{2}-h]$, and $\mathbf{v}$ ranges over all grid points of a grid of $S^{d-1}$ and grid width converging to zero with increasing sample size.

The testing procedure is as follows. For any $\mathbf{b}_{0}$ in $(\mathbf{a}_{1},\mathbf{a}_{2}),$ let $\mathcal{T}_{n}^{\mathbf{b}_{0}}$ be the set of all sequences of triples $(\mathbf{t}_{n},h_{n},\mathbf{v}_{n})_{n}\in\mathcal{T}_{n}$ such that $ch_{n}\geq\|\mathbf{b}_{0}-\mathbf{t}_{n}\|\geq 2h_{n}$ for some sufficiently large $c>2$ and $\angle(\mathbf{t}_{n}-\mathbf{b}_{0},\mathbf{v}_{n})\rightarrow 0$ as $n\rightarrow\infty$, where $\angle$ denotes the angle between two vectors. The previous conditions ensure that several such sequences can be found. If for all triples in $\mathcal{T}_{n}^{\mathbf{b}_{0}}$ all local tests (3.16) reject the hypotheses (3.5), we have evidence for the existence of a mode at the point $\mathbf{b}_{0}$. By choosing the test locations as the vertices of an equidistant grid no prior knowledge about the location of $\mathbf{b}_{0}$ has to be assumed. Theorem 3.4 below states that the procedure detects all modes of the density with probability converging to one as $n\rightarrow\infty$.

Theorem 3.4**.**

Assume the conditions of Theorem 3.2 and suppose that for any mode $\mathbf{b}_{0}$ in $(\mathbf{a}_{1},\mathbf{a}_{2})$ there are functions $g_{\mathbf{b}_{0}}:\mathbb{R}^{d}\rightarrow\mathbb{R}$, $r_{\mathbf{b}_{0}}:\mathbb{R}\rightarrow\mathbb{R}$ such that the density $f_{\bm{\beta}}$ has a representation of the form

[TABLE]

in an open neighborhood of $\mathbf{b}_{0}$. Furthermore, let $g_{\mathbf{b}_{0}}$ be differentiable in an open neighborhood of $\mathbf{b}_{0}$ with $g_{\mathbf{b}_{0}}(\mathbf{b})=o(1)$ and $\langle\nabla g_{\mathbf{b}_{0}}(\mathbf{b}),\mathbf{e}\rangle=o(\|\mathbf{b}-\mathbf{b}_{0}\|)$ when $\mathbf{b}\rightarrow\mathbf{b}_{0}$ for all $\mathbf{e}\in\mathbb{R}^{d}$ with $\|\mathbf{e}\|=1$. In addition, let $\widetilde{f}_{\mathbf{b}_{0}}$ be differentiable in an open neighborhood of zero with $\widetilde{f}_{\mathbf{b}_{0}}^{\prime}(h)\leq-ch(1+o(1))$ for $h\rightarrow 0$.

If $\min\{h:h\in H(\mathcal{T}_{n})\}\geq C\log(n)^{1/(2d+3)}n^{-1/(2d+3)}$ is nonempty for some sufficiently large constant $C>0,$ then the procedure described in the previous paragraph detects the mode $\mathbf{b}_{0}$ with probability converging to one as $n\rightarrow\infty.$

The rate for the localization of the modes is $n^{-1/(2d+3)}$ (up to some logarithmic factor). Since $2d+3=d+2\tfrac{d-1}{2}+4,$ this rate matches the optimal rate for mode detection in an inverse problem with ill-posedness of degree $(d-1)/2$ over a 2-Hölder class. Assumption 4 requires, however, $h_{\,\textnormal{max}}\lesssim\log(n)^{-14\gamma/(d-1)-5}n^{2\gamma/(d-1)-1}.$ To be able to include scales of the order $n^{-1/(2d+3)}$ we need $\gamma>(d^{2}-1)/(2d+3).$ The right hand side is smaller than one for $d=2,3.$

4 Finite sample properties

In this section we illustrate the finite sample properties of the proposed test in a bivariate and a trivariate setting. In the bivariate setting we illustrate how simultaneous tests for the hypotheses (3.3) and (3.4) can be used to obtain a graphical representation of the local monotonicity properties of the density. In the trivariate setting we investigate the performance of the test for modality at a given point $\mathbf{b}_{0}$ (see the hypotheses in (3.5)) and the dependence of its power on the distribution of $\mathbf{X}$.

As test function we consider the simplest polynomial which satisfies the conditions of Assumption 1 for $d=2,3$, that is,

[TABLE]

with $c$ such that $\int\phi=1.$ Figure 2 displays the function $\mathcal{H}_{d}\widetilde{\phi}^{(d-1)}$ for $d=2,3.$ Throughout this section the nominal level is fixed as $\alpha=0.05$, and all level and power statements are in percent. Except for Table 2, none of the simulations in this section assumes knowledge of the design density $f_{\bm{\Theta}}$ or uses a parametric specification of it.

4.1 Inference about local monotonicity of a bivariate density

We follow the multiscale approach in Section 3.1 to obtain a graphical representation of the monotonicity behavior for a bivariate density of random coefficients. To test the hypotheses (3.4) we use (3.16) with $\mathcal{T}_{n}=\mathcal{T}_{\mathbf{t}}\times\{h_{0}\}\times\mathcal{T}_{\mathbf{v}}.$ Here, $h=h_{0}=0.5$ is fixed and the set of test locations $\mathcal{T}_{\mathbf{t}}$ is defined as the set of vertices on an equidistant grid in the square $[-1,2]^{2}$ with width one. Finally, the set of test directions is

[TABLE]

The data are simulated with $f_{\bm{\beta}}$ the density of the normal mixture $\tfrac{1}{3}{\cal N}((-0.4,-0.57)^{\top},0.2I)+\tfrac{1}{3}{\cal N}((1.5,-0.52)^{\top},0.2I)+\tfrac{1}{3}{\cal N}((0.45,1.6)^{\top},0.15I).$ The design is chosen such that ${\bm{\Theta}}_{i}$ is uniformly distributed on the sphere $\mathbb{S}^{1}$. Figure 1 in Section 3.1 displays the monotonicity behavior of the density $f_{\bm{\beta}}$ based on sample size $n=20000$. Each arrow at a location $\mathbf{t}$ in direction $\mathbf{v}$ displays a rejection of a hypothesis (3.4). The map indicates the existence of modes around the points $(-0.5,-0.5)^{\top}$, $(1.5,-0.5)^{\top}$, and $(0.5,1.5)^{\top}$ and thus detects the true modes fairly well.

4.2 Influence of the design on the power

Given the random coefficients model with $d=3$, we study the power of the test for the existence of a mode at a given location $\mathbf{b}_{0}$ considering only few local tests. The postulated mode is given by the point $\mathbf{b}_{0}=(0,0,0)^{\top}$ and we take $\mathcal{T}_{n}=\{(\mathbf{t},h,\mathbf{v}):h=1,\mathbf{t}=\mathbf{v}\in\{\pm\mathbf{e}_{i},i=1,2,3\}\}$ with $\mathbf{e}_{i}$ the $i$-th standard unit vector in $\mathbb{R}^{3}.$ We conclude that $f_{\bm{\beta}}$ has a local maximum at the point $\mathbf{b}_{0}$, whenever all hypotheses $H_{0,-}^{\mathbf{t}_{j},h_{0},\mathbf{v}_{j}},$ $j=1,\ldots,6,$ are rejected, i.e.

[TABLE]

where $\kappa^{\mathbf{t}_{j},h_{0},\mathbf{v}_{j}}_{n}(\alpha)$ is defined by (3.14). Recall that the quantiles $\kappa^{\mathbf{t}_{j},h_{0},\mathbf{v}_{j}}_{n}(\alpha)$ are constructed in such a way that the probability of at least one false rejection within the six tests (4.1) is bounded by $\alpha$. However, the mode test detects the presence of a mode whenever all six tests (4.1) are rejected at the same time. The multiscale method is therefore rather conservative for the specific task of mode detection. In this simulation, we also study a calibrated version of the test where the quantiles are chosen such that the test keeps its nominal level $\alpha=0.05$ and detects the presence of a non existing mode in about 5 percent of the simulation runs. For the calibration of the test we work under the null hypothesis assuming that $f_{{\bm{\beta}}}$ is uniform. Therefore, knowledge about the true unknown density $f_{\bm{\beta}}$ is not required.

Numerical simulations for random coefficients model without intercept: At first, we consider model (1.1) with uniform design $\mathbf{X}_{i}\sim\operatorname{Unif}[-5,5]^{3}.$ To study the level of the test, we used $f_{\bm{\beta}}({\bm{\beta}})\propto\mathbf{1}({\bm{\beta}}\in[-5,5]^{3}).$ For the power we took $f_{\bm{\beta}}$ as the density of a trivariate standard normal distribution. All results are based on the local tests (4.1) and 1000 repetitions. Level simulations with the theoretical quantiles confirm that the multiscale test keeps its nominal level as the percentage of false rejections of at least one of the six hypotheses in (4.1) for sample size $n\in\{250,500,1000\}$ is nearly $5$ percent. The results of mode test are reported in Table 1.

Next, we investigate an asymmetric distribution of the directions ${\bm{\Theta}}_{i}$ by sampling $\mathbf{X}_{i}\sim\mathcal{N}((3,0,0)^{\top},2I)$, with $I$ the $3\times 3$ identity matrix. We only consider the calibrated mode test. Results are reported in columns five and six in Table 1. Compared to the case of uniform design, we observe a decrease in the power of the mode test (4.1). The explanation is that the uniform design on the $\mathbf{X}_{i}$ induces a more uniform distribution of ${\bm{\Theta}}_{i}$ on the sphere which makes it simpler to recover information about the joint density as discussed in Section 2.

Numerical simulations for random coefficients model with intercept: We study model (1.2) with $d=3.$ In a first simulation, we sample the random vectors $(X_{i,2},X_{i,3})^{\top}$ from a standard bivariate Cauchy distribution, such that the density $f_{\bm{\Theta}}$ is constant. Except for the different design, we consider otherwise the same test settings as above. The simulated level and power of the calibrated version of the test (4.1) are reported in Table 2. To investigate the influence of the estimation of the design density $f_{\bm{\Theta}}$ on the power of the test we also perform simulations in which we assume that the density $f_{\bm{\Theta}}$ is known to be constant. These are shown in Table 2, fourth and fifth column.

Compared to the power approximations for unknown $f_{\bm{\Theta}}$, we observe only a slight increase in power of the test for known $f_{\bm{\Theta}}$.

Finally, we consider two designs which do not satisfy Assumption 3. Table 3 reports the level and power for the same setting as above except that now $(X_{i,2},X_{i,3})^{\top}$ is drawn from a standard normal distribution or $(X_{i,2},X_{i,3})^{\top}\sim\operatorname{Unif}[-5,5]^{2}.$ We observe only a slight decrease in the power of the test for normally distributed design compared to the setting where Assumption 3 holds. Even under uniform design, the test performs fairly well.

4.3 Multiscale mode testing

For multimodal densities which have a second mode close to the test location $\mathbf{b}_{0}$ testing different bandwidths simultaneously can be advantageous to separate the modes. This is illustrated by the following example, where we consider the random coefficients model without intercept with $d=2.$ The data are simulated with $f_{\bm{\beta}}$ being the density of the normal mixture

[TABLE]

We consider a design such that ${\bm{\Theta}}$ is uniformly distributed on the circle $\mathbb{S}^{1}$. The level $\alpha$ is fixed to five percent. We conducted simultaneously twelve tests with three different scales $h\in\{0.5,1,2.5\}$ for the hypotheses (3.5) with $\mathbf{b}_{0}=(0,0)^{\top}$. The tests are given by $\{(\mathbf{t}_{i},h_{j},\mathbf{v}_{i}):i=1,\ldots,4;h_{j}\in\{0.5,1,2.5\},\mathbf{v}_{i}\in\{\pm\mathbf{e}_{1},\pm\mathbf{e}_{2}\},\mathbf{t}_{i}=h_{j}\mathbf{v}_{i}\}.$ We analyze the outcome of the twelve tests in two ways. Firstly, we investigate the performance of each of the three tests for modality for the bandwidths $h=2.5$, $h=1$ and $h=0.5$ separately. Secondly, we consider the performance of a combined test for modality using two scales $h=1$ and $h=0.5$. The power approximations for sample sizes $n\in\{2000,5000,15000\}$ based on 1000 repetitions are displayed in Table 4.

Figure 3 illustrates the results of the twelve tests for the hypotheses (3.5) conducted simultaneously. Each arrow at a location $\mathbf{t}$ in direction $\mathbf{v}$ displays a rejection of a hypothesis (3.5) and the length of the arrows corresponds to the respective bandwidths.

The results in Table 4 show that the mode test for $h=2.5$ detects in all cases a mode in a neighborhood of $(0,0)^{\top}$. However, the bandwidth $h=2.5$ is too large to distinguish between the underlying modes of $f_{\bm{\beta}}$ at $(0,0)^{\top}$ and $(2,0)^{\top}$. The test for the bandwidth $h_{0}=1$ fails to detect the mode as $h=1$ is still too large to separate the two modes of $f_{\bm{\beta}}.$ On the contrary, $h_{0}=0.5$ is too small to detect the decrease with small slope corresponding to the eigenvalue $0.4$ of the covariance matrix in the first mixture component of $f_{\bm{\beta}}.$ Note that this effect vanishes for increasing sample size. By conducting the tests for the bandwidths $h=1$ and $h=0.5$ simultaneously, we are able to detect the mode at $(0,0)^{\top}$ in most of the simulations.

4.4 Comparison with the parametric random coefficient model

It is common in the literature to use a parametric specification of the random coefficients. Usually, ${\bm{\beta}}$ is assumed to have a multivariate normal distribution. In this section, we compare our nonparametric approach with mode estimation under a parametric specification of the random coefficient distribution. Our findings are similar to what is usually observed when parametric and nonparametric methods are compared. If the parametric specification of the model is consistent with the data generating process, the parametric approach outperforms the nonparametric one in terms of estimation errors, power of tests, and computation time. However, when the parametric specification and the data generating process differ, a model misspecification bias is present in the parametric estimation. This will not vanish, even if the sample size is large. The nonparametric model does not suffer from this problem and can in these cases perform better than the parametric model. We illustrate this effect in the following simulations.

We consider the parametric random coefficients model $Y=\tilde{\bm{\beta}}^{\top}\mathbf{X}$ where the density of $\tilde{\bm{\beta}}$ belongs to some parametric family $f_{\tilde{\bm{\beta}}}(\mathbf{b};\eta)$ with parameter $\eta$. Note that we can rewrite the equation as a mixed model $Y={\bm{\gamma}}^{\top}\mathbf{X}+\mathbf{u}^{\top}\mathbf{X}$. Here ${\bm{\gamma}}=\mathbb{E}(\tilde{\bm{\beta}})$ are fixed effects and $\mathbf{u}=\tilde{\bm{\beta}}-\mathbb{E}(\tilde{\bm{\beta}})$ are random effects with expectation $\bm{0}$. If $f_{\tilde{\bm{\beta}}}(\mathbf{b};\eta)$ is unimodal with a mode at $\mathbb{E}(\tilde{\bm{\beta}})$, we just need to derive confidence statements for ${\bm{\gamma}}$. This is in particular true if $f_{\tilde{\bm{\beta}}}(\mathbf{b};\eta)$ is a normal density, which is the most common choice. A simple and efficient way to estimate ${\bm{\gamma}}$ is heteroscedasticity robust linear regression. Obviously, the computational complexity of this method is much smaller than the complexity of our algorithm.

Our test example is the random coefficients model with intercept (1.2) with $d=3$ and $(X_{i,2},X_{i,3})^{\top}$ sampled from a standard bivariate Cauchy distribution. The distribution of the random coefficients ${\bm{\beta}}_{i}$ in the data generating process is given by (i) a standard Gaussian, (ii) the normal mixture $0.5\mathcal{N}(\bm{0},0.1I)+0.5\mathcal{N}((2,0,0)^{\top},0.1I)$ and (iii) an exponential-2 distribution for $\beta_{i,1}$ and $(\beta_{i,2},\beta_{i,3})^{\top}\sim\mathcal{N}(\bm{0},0.1I)$ independent of the first component. Table 5 reports estimates for ${\bm{\gamma}}$ obtained by transforming $Y$ and $\mathbf{X}$ to $S$ and ${\bm{\Theta}}$ and running a heteroscedasticity robust regression.

In column (i) of Table 5 the parametric assumption holds and the procedure detects the true mode $\bm{0}$ of the density with high precision compared to the bandwidth choice $h=1$ in the simulations presented in Table 2 in Section 4.2. In contrast, for the bimodal density in column (ii) the coefficient vector $(1,0,0)^{\top}$ does not describe a representative member of the population because the misspecification bias is too large. For the skewed distribution of column (iii) the estimator also fails to detect the mode of the density for the same reason.

By applying our testing procedure in the setting of columns (ii) and (iii) we can show that the OLS results do not represent modes of the density. To this end, we set in (ii)

[TABLE]

and in (iii)

[TABLE]

Both tests reject $H_{0,-}^{\mathbf{t},h,\mathbf{v}}$ (in (ii) in $100\%$ and in (iii) in $95.5\%$ percent of 1000 repetitions). Therefore, our procedure shows that neither $(1,0,0)^{\top}$ in (ii) nor $(2,0,0)^{\top}$ in (iii) are modes of the underlying density. Of course, the parametric model would be able to detect the modes in cases (ii) and (iii) with a different parametric specification. However, this would require considerable a priori knowledge about the data. If we did not interpret the results as estimators for the mode but as estimators for $\mathbb{E}({\bm{\beta}})$, the parametric method would perform well.

5 Application to consumer demand data

Heterogeneity of consumers is a major challenge in modeling and estimating consumer demand. In several different demand models random coefficients were proposed to account for the heterogeneity in the population of consumers.

5.1 Model and data

In this section we are interested in the almost ideal demand system (AIDS) which was initially proposed by Deaton and Muellbauer, (1980) with fixed coefficients. This model does not explain demand for a product itself but explains the budget share spent on a product by a linear equation. The explanatory variables are log prices and the log of total expenditure divided by a price index. A detailed discussion of the model is contained in Lewbel, (1997).

Fixed coefficients in this model mean that all consumers are assumed to react in the same way when the price of a product changes. It is well known that some consumers are very price sensitive and change their behavior significantly with small variations in prices while other consumers are less price sensitive. This type of heterogeneity can be modeled by a random coefficient on log prices which is assumed to vary across the population of consumers. A similar argument suggests a random coefficient on log total expenditure. Recently, applications of the AIDS using a nonparametric random coefficient specification instead of fixed coefficients were presented in Hoderlein et al., (2015) and Breunig and Hoderlein, (2018).

We apply our multiscale test to detect modes in the random coefficients for budget shares for food at home (BSF)

[TABLE]

Food expenditure is a large fraction of total expenditure and is roughly about 20%.

We analyze the data of the British Family Expenditure Survey which ran from 1961 to 2001. It reported yearly cross sections for household income, expenditure and other characteristics of roughly 7000 households. We use data of the years 1997–2001 only which gives a sample size of about 33000. Budget shares of food are generated by dividing the expenditure for all food by total expenditure. Food prices are reported as relative prices in comparison to a general prize index. The variable $TotExp$ is normalized to January 2000 real prizes.

Assumption 3 and the numerical simulations in Section 4 suggest that our test has more power when the normalized regressors are approximately uniform on the sphere. We can achieve this by symmetrizing the design in model (4.1) as follows:

[TABLE]

The relation of the modified model to the random coefficients in (4.1) is $\beta_{i,1}=\tilde{\beta}_{i,1}-5\tilde{\beta}_{i,2}-0.3\tilde{\beta}_{i,3}$, $\beta_{i,2}=\tilde{\beta}_{i,2}$, $\beta_{i,3}=25\tilde{\beta}_{i,3}$. Observations of the new variable $\ln(TotExp_{i})-5$ lie between $-5$ and $3.7$. The observations of $25\ln(FoodPrice_{i})-0.3$ range from $-1$ to $1.3$.

5.2 Results

For a first evaluation of the data we assumed fixed coefficients in model (4.2) and estimated the model with ordinary least squares (OLS).

In order to find modes of the density we conducted simultaneously tests on the $5\%$ level of the form (3.5) on the two scales $h_{1}=0.75$ and $h_{2}=0.5$. Recall from Section 4.2 that our testing procedure also performs well when the $\mathbf{X}_{i}$ are not Cauchy distributed. Therefore, to obtain a testing procedure which is more flexible with respect to the design, we use the nonparametric density estimator $\widehat{f}_{\bm{\Theta}}$ instead of a parametric estimation procedure. We were testing for modes on the equidistant grid covering $[-1,1]^{3}$ with grid width $1$. Hence, the grid had 27 nodes. For every grid point $\mathbf{b}\in\mathbb{R}^{3}$ tests of the hypotheses (3.5) were conducted for the directions and locations

[TABLE]

where $\mathbf{e}_{1},\mathbf{e}_{2},\mathbf{e}_{3}\in\mathbb{R}^{3}$ denote the standard unit vectors of $\mathbb{R}^{3}$. We detected a single mode in the neighborhood of the grid point $(0,0,0)^{\top}$ for the tests with bandwidth $h_{1}=0.75$. The test for the bandwidth $h_{2}=0.5$ did not detect a mode.

In the following we use nonparametric density estimation to motivate hypotheses for the testing procedure. It is important that this estimate and the test are independent, otherwise the testing procedure would be biased and could not guarantee a bound on the error rate. We meet the requirement by splitting the sample in two independent equally sized sub-samples. The first sub-sample is used for nonparametric estimation of the random coefficient density in model (4.2) with the estimator in Hoderlein et al., (2010). Figure 4 gives contour plots for the joint densities of $f_{\tilde{\beta}_{1},\tilde{\beta}_{2}}$, $f_{\tilde{\beta}_{1},\tilde{\beta}_{3}}$, and $f_{\tilde{\beta}_{2},\tilde{\beta}_{3}}$ based on about 16500 observations. We chose the smoothing parameters $h$ and $g$ in the estimator in Hoderlein et al., (2010) equal to 0.05 and 0.1, respectively. Note that these bandwidth choices do not affect the level of the test performed below as the test does not depend on this estimator. The nonparametric estimate suggests that the random coefficient density of $f_{\tilde{\beta}_{1},\tilde{\beta}_{2},\tilde{\beta}_{3}}$ has one (well-pronounced) mode close to

[TABLE]

This is consistent with the results of the test above which found a mode close to $(0,0,0)$. Since the marginal densities of $\beta_{1},\beta_{2},\beta_{3}$ are nearly symmetric it is also consistent that the mode is close to the OLS estimates given in Table 6. With a significantly skewed or with a multimodal random coefficient density location of modes would differ from OLS.

With the second sub-sample we studied whether a mode can be found for the location given by (4.3) on a smaller scale than in the test above. This would then indicate that the location of the mode is not $(0,0,0)^{\top}$. Our primary interest is in the coefficients on $\tilde{\beta}_{2}$ and $\tilde{\beta}_{3}$ on total expenditure and food prices. In order to see if the mode is indeed in a location where $\tilde{\beta}_{2}<0$ and $\tilde{\beta}_{3}>0$, we conduct two mode tests (4.3) simultaneously for the bandwidths $h_{1}=0.07$ and $h_{2}=0.02$, that is, we test twelve hypotheses (3.5) for the following locations and directions $(\mathbf{t}_{i},\mathbf{v}_{i}).$ On scale $h_{1}=0.07$, we tested

[TABLE]

and on scale $h_{2}=0.02,$

[TABLE]

The test rejected all local hypotheses on scale $h_{1}=0.07$ but not on scale $h_{2}=0.02.$ This gives evidence that the mode is in a location where $\tilde{\beta}_{2}$ is negative but we cannot decide whether $\tilde{\beta}_{3}$ is positive at the mode.

Let us return to the initial model (4.1). The results of our test give evidence that a mode exists close to

[TABLE]

with strong evidence that $\beta_{2}$ is indeed negative. This vector of coefficients describes a representative member of the majority of consumers. It suggests that in the majority group food budget shares decrease with increasing log total expenditure. The nonparametric estimate in Figure 4 shows that there is considerable variance among consumers around this representative member.

Acknowledgments

J. Schmidt-Hieber was partially funded by a TOP II grant from the Dutch science foundation. F. Dunker acknowledges support by the Ministry of Education and Cultural Affairs of Lower Saxony in the project Reducing Poverty Risk. K. Proksch acknowledges financial support by the German Research Foundation DFG through subproject A07 of CRC 755. K. Eckle has been supported by the Collaborative Research Center “Statistical modeling of nonlinear dynamic processes” (SFB 823, Project C4) of the German Research Foundation (DFG). We are very grateful to two referees and an associate editor for their constructive comments, which led to substantial improvement of an earlier version of this manuscript.

Supplementary material to “Tests for qualitative features in the random coefficients model”

Appendix A Nonparametric estimators for the densities $f_{\mathbf{\Theta}}$ and $f_{S,\mathbf{\Theta}}$

In this section we discuss the estimation of the densities $f_{\bm{\Theta}}$ and $f_{S,{\bm{\Theta}}}$ and related quantities. We use kernel density estimators based on the second half of the observations $(S_{i},{\bm{\Theta}}_{i}),$ $i=n+1,\ldots,2n$. The density of the random vector $\mathbf{X}_{1}$ is denoted by $f_{\mathbf{X}}$. In the random coefficients model without intercept $f_{\mathbf{X}}$ is a $d$-variate density and a $(d-1)$-variate density in the random coefficients model with intercept. Throughout the following, $K:\mathbb{R}\rightarrow\mathbb{R}$ is assumed to be Lipschitz continuous, non-negative and $\int K=1.$

In the random coefficients model without intercept, we introduce the kernel density estimator

[TABLE]

with normalization constant

[TABLE]

As shown in Bai et al., (1988), the integral does not depend on ${\bm{\theta}}$ and $C(h_{*})$ converges to some positive constant as $h_{*}\rightarrow 0$. For the joint density of $(S,{\bm{\Theta}}),$ we propose the kernel density estimator

[TABLE]

In the random coefficients model with intercept the symmetrizations $S_{i}=\zeta_{i}Y_{i}/\|\mathbf{X}_{i}\|$ and ${\bm{\Theta}}_{i}=\zeta_{i}\mathbf{X}_{i}/\|\mathbf{X}_{i}\|$ with Rademacher variables $\zeta_{i}$ correspond to point reflections of the densities at the origin. Thus, the density $f_{\bm{\Theta}}$ is in general not continuous on the boundary of the hemisphere $\mathbb{S}^{d-1}_{+}$ (see also (3.11)). Smoothness is, however, necessary to control the bias. Therefore, we use a two step procedure for the estimation of $f_{\bm{\Theta}}$ and $f_{S,{\bm{\Theta}}}$ in the random coefficients model with intercept. First, we estimate the density of the non-symmetrized samples $Y_{i}/\|\mathbf{X}_{i}\|$ and $\mathbf{X}_{i}/\|\mathbf{X}_{i}\|$, $i=n+1,\ldots,2n$ on the hemisphere $\mathbb{S}^{d-1}_{+}$ and on $\mathbb{R}\times\mathbb{S}^{d-1}_{+}$, respectively. The estimators for $f_{\bm{\Theta}}$ and $f_{S,{\bm{\Theta}}}$ are then the same as (4.1) and (4.2) except that now ${\bm{\theta}}\in\mathbb{S}^{d-1}_{+}$ and the normalization constant $C(h_{*})$ is replaced by a function ${\bm{\theta}}\mapsto C(h_{*},{\bm{\theta}})$, defined by

[TABLE]

In a second step, we symmetrize the estimators and divide by two to get estimators of the densities $f_{\bm{\Theta}}$ and $f_{S,{\bm{\Theta}}}$ on the whole domain.

In the next lemma we establish convergence of these estimators. If in Lemma A.1 $f_{\mathbf{X}}$ follows a multivariate Cauchy distribution, then $\gamma:=1$. Otherwise, $\gamma$ comes from Assumption 3 (iii).

Lemma A.1**.**

Suppose Assumption 3 is satisfied for some $\gamma>0$ and set $\gamma=1$ in the case of the random coefficients model with intercept. In both models, the estimator $\widehat{f}_{\bm{\Theta}}$ with bandwidth $h_{*}=O(\log(n)^{-3/\gamma})$ and $h_{*}\geq\log(n)^{7/(d-1)}n^{-1/(d-1)}$ satisfies

[TABLE]

The proof is delayed to the end of this section. Let us now discuss properties of the density $f_{S,{\bm{\Theta}}}(s,{\bm{\theta}}).$ By (2.1), $f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})=f_{\bm{\Theta}}({\bm{\theta}})Rf_{\bm{\beta}}(s,{\bm{\theta}}).$ Under Assumption 2, $f_{\bm{\beta}}$ is compactly supported and consequently, $f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})$ and $s\mapsto\log(|s|)^{2}f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})$ for $|s|\geq 1$ are uniformly bounded. Moreover, $f_{S,{\bm{\Theta}}}$ is Hölder continuous with Hölder constant $\gamma$. This is a straightforward consequence of the Hölder $\gamma$-continuity of $f_{\bm{\Theta}}$ shown in the proof of Lemma A.1 and the identity

[TABLE]

where ${\bm{\theta}}_{1},\ldots,{\bm{\theta}}_{d-1}$ denote an orthonormal basis of the orthogonal complement of $\rm{span}\{{\bm{\theta}}\}$, together with the compact support and the Lipschitz-continuity of $f_{\bm{\beta}}$ (following from Assumption 2). Moreover, the properties of $f_{S,{\bm{\Theta}}}$ discussed above also hold in quantum homodyne tomography under Assumption 2’. We point out that the marginal densities of the Wigner function, which are given by the Radon transform, are nonnegative. This will be used later to bound the standard deviation $\sigma_{\mathbf{t},h,\mathbf{v}}$ away from zero.

If in Lemma A.2 $f_{\mathbf{X}}$ follows a multivariate Cauchy distribution, then $\gamma:=1$. Otherwise, $\gamma$ comes from Assumption 3 (iii).

Lemma A.2**.**

Let Assumptions 2 resp. 2’ and 3 hold and consider the estimator $\widehat{f}_{S,{\bm{\Theta}}}$ in (4.2) with bandwidth choice $h_{+}=O(\log(n)^{-3/\gamma})$ and $h_{+}\geq\log(n)^{3/d}n^{-1/(2d)}$. Then

[TABLE]

For the estimation of the test statistic $\widehat{T}_{\mathbf{t},h,\mathbf{v}}$ and the limiting process $\widehat{X}_{\mathbf{t},h,\mathbf{v}}$ the quantities $1/f_{\bm{\Theta}}$ and $\sqrt{f_{S,{\bm{\Theta}}}}$ need to be estimated. The functions $(\cdot)^{-1}$ and $\sqrt{\cdot}$ are not smooth in zero and we therefore introduce the cut-off estimators

[TABLE]

By the boundedness from below of $f_{\bm{\Theta}}$ and Lemma A.1 it holds $\widetilde{f}_{\bm{\Theta}}=\widehat{f}_{\bm{\Theta}}$ almost surely for $n$ sufficiently large.

Proof of Lemma A.1.

We only consider the case without intercept, that is, (i)-(iii) in Assumption 3 hold. In the case with intercept, we use $C(h_{*})\leq C(h_{*},{\bm{\theta}})\leq 2C(h_{*})$ and the fact that $C(h_{*},{\bm{\theta}})$ is Lipschitz-continuous with respect to ${\bm{\theta}}$ to arrive at the same conclusion.

To prove (i) observe that

[TABLE]

Here, we used the compact support of $K$ and the identity $1-\langle{\bm{\theta}}^{\prime},{\bm{\theta}}\rangle=\|{\bm{\theta}}^{\prime}-{\bm{\theta}}\|^{2}/2$. Since $f_{\bm{\Theta}}({\bm{\theta}})=\int_{0}^{\infty}r^{d-1}f_{\mathbf{X}}(r{\bm{\theta}})dr$, we have by Assumption 3 (iii) $|f_{\bm{\Theta}}({\bm{\theta}}^{\prime})-f_{\bm{\Theta}}({\bm{\theta}})|\leq\int_{0}^{\infty}r^{d-1}|f_{\mathbf{X}}(r{\bm{\theta}}^{\prime})-f_{\mathbf{X}}(r{\bm{\theta}})|dr\lesssim h_{*}^{\gamma}$ for $\|{\bm{\theta}}^{\prime}-{\bm{\theta}}\|\lesssim h_{*}$. By definition of the constant $C(h_{*})$, we obtain \big{|}\mathbb{E}\big{[}\widehat{f}_{\bm{\Theta}}({\bm{\theta}})\big{]}-f_{\bm{\Theta}}({\bm{\theta}})\big{|}\lesssim h_{*}^{\gamma} and this proves (i).

Next, we bound the stochastic error term (ii) using an entropy argument and Bernstein’s inequality. Observe that by the Lipschitz-continuity of $K,$ $\widehat{f}_{\bm{\Theta}}-\mathbb{E}[\widehat{f}_{\bm{\Theta}}]$ is Lipschitz-continuous with Lipschitz constant of order $h_{*}^{-d-1}.$ For $a_{n}:=\sqrt{\frac{\log(n)}{nh_{*}^{d-1}}},$ let $\{{\bm{\theta}}_{j}:j=1,\ldots,M\}$ be defined as the set of smallest cardinality such that $\bigcup_{j=1}^{M}B_{c^{\prime}h_{*}^{d+1}a_{n}}({\bm{\theta}}_{j})\supset\mathbb{S}^{d-1}$ for some constant $c^{\prime}>0.$ If $c^{\prime}>0$ is chosen small enough, then

[TABLE]

In order to bound the probability, we apply Bernstein’s inequality to $\widehat{f}_{\bm{\Theta}}({\bm{\theta}}_{j})-\mathbb{E}[\widehat{f}_{\bm{\Theta}}({\bm{\theta}}_{j})]=\sum_{i=n+1}^{2n}Z_{i}$ with

[TABLE]

We find $|Z_{i}|\leq\frac{C_{1}}{nh_{*}^{d-1}}$ for some constant $C_{1}>0,$ and for some constant $C_{2}>0,$

[TABLE]

using the boundedness of $K$ and $f_{\bm{\Theta}}$ as well as the definition of $C(h_{*})$. Hence, an application of Bernstein’s inequality yields with (4.4),

[TABLE]

Since $M$ is a polynomial power of $n,$ the claim follows by choosing the constant $c$ large enough.

For (iii) one proceeds similarly with the choice $a_{n}=(\log(n)\log\log(n))^{-1}$ using the summability of the probabilities. ∎

To prove Lemma A.2, we make use of a slightly modified version of Proposition 2.2 in Giné and Guillou, (2001) which we state below as Proposition A.1. If $\mathcal{F}$ is a uniformly bounded class of measurable functions on a measurable space $(S,\mathcal{S})$ with a measurable and bounded envelope $F$, then $\mathcal{F}$ is said to be a measurable uniformly bounded VC class of functions if there are constants $A,v>0$ such that

[TABLE]

for all $0<\varepsilon<1$, where $N(T,d,\varepsilon)$ denotes the $\varepsilon$-covering number of the metric space $(T,d)$ and the supremum is taken over all probability measures on $(S,\mathcal{S})$.

Proposition A.1**.**

Let $P$ be any probability measure on $(S,\mathcal{S})$ and let $\xi_{i},\;i=1,\ldots,n,$ be independent with common law $P$. Let further $\mathcal{F}$ be a measurable uniformly bounded VC class of functions and let $\sigma^{2}$ and $U$ be any numbers such that $\sigma^{2}\geq\sup_{f\in\mathcal{F}}\operatorname{Var}_{P}f,\;U\geq\sup_{f\in\mathcal{F}}\|f\|_{\infty}$ and $0<\sigma\leq U$. Then there exist universal constants $C,K^{\prime},L>0$ such that the exponential inequality

[TABLE]

is valid for all $t\geq C(vU\log(AU\sigma^{-1})+\sqrt{vn}\sigma\sqrt{\log(AU\sigma^{-1})}).$

In contrast to Proposition 2.2 in Giné and Guillou, (2001), Proposition A.1 contains the explicit dependence of the right hand side of (LABEL:1.1) on the constants $A$ and $v$.

Proof of Lemma A.2.

Similarly as in the proof of Lemma A.1, it is enough to consider the random coefficients model without intercept only and to work under Assumption 3, (i)-(iii). If the design density is multivariate Cauchy, we can derive the properties in a similar way for $\mathbb{R}\times\mathbb{S}^{d-1}_{+}$. An upper bound of the bias can be derived similarly to Lemma A.1 (i). For the stochastic error

[TABLE]

we apply Proposition A.1 to the function class

[TABLE]

which depends on $n$ via $h_{+}$. By the boundedness of $K$ we find $U,\sigma\lesssim 1$. To show that $\mathcal{F}_{n}$ is a VC class of functions, we introduce a discretization of $\mathbb{R}\times\mathbb{S}^{d-1}$ as follows: Let $c>0$ be a sufficiently small constant only depending on the kernel $K$. We chose a grid $\{{\bm{\theta}}_{j}:j=1,\ldots,M_{1}\}$ of $\mathbb{S}^{d-1}$ with grid width at most $c\varepsilon h_{+}^{2}$. Obviously, this is possible with $M_{1}\lesssim(\varepsilon^{-1}h_{+}^{-2})^{d-1}$. Moreover, introduce the set of intervals $I_{k}=[k,k+1)$, $k\in\mathbb{Z}$. For each probability measure $Q$ there are at most $\lceil(c\varepsilon)^{-2}\rceil$ sets $I_{i_{j}}\times\mathbb{S}^{d-1},\;j=1,\ldots,\lceil(c\varepsilon)^{-2}\rceil,$ such that $Q(I_{i_{j}}\times\mathbb{S}^{d-1})\geq(c\varepsilon)^{2}$. Let $\{s_{j}:j=1,\ldots,M_{2}\}$ be an equidistant grid of

[TABLE]

with grid width $ch_{+}\varepsilon$ and let $s_{M_{2}+1}$ denote an arbitrary point in $\widetilde{I}_{h_{+}}^{C}.$ Basic calculations show $M_{2}\lesssim\varepsilon^{-3}h_{+}^{-1}.$ Moreover, the subset of $\mathcal{F}_{n}$ indexed by

[TABLE]

is an $\varepsilon$-covering set of $\mathcal{F}_{n}$. To see this, fix $(s,{\bm{\theta}})\in\widetilde{I}_{h_{+}}\times\mathbb{S}^{d-1}$. Then

[TABLE]

by the Lipschitz continuity of $K$. Hence, by construction of the set $\widetilde{I}_{h_{+}}\times\mathbb{S}^{d-1}$ there exists $j\in\{1,\ldots,M_{1}(M_{2}+1)\}$ such that

[TABLE]

For $(s,{\bm{\theta}})\in(\widetilde{I}_{h_{+}}\times\mathbb{S}^{d-1})^{C}$ we obtain

[TABLE]

since the support of $K(\frac{\cdot-s}{h_{+}})$ is compact and does not intersect with any of the sets $I_{i_{j}}\times\mathbb{S}^{d-1},\;j=1,\ldots,\lceil(c\varepsilon)^{-2}\rceil$ for $h_{+}$ sufficiently small. A similar argument applies to $K(h_{+}^{-1}(\cdot-s_{M_{2}+1}))$. Hence,

[TABLE]

and $\mathcal{F}_{n}$ is a VC class of functions with $v=d+2$ and $A=A_{n}=h_{+}^{(-2d+1)/(d+2)}$. An application of Proposition A.1 yields

[TABLE]

for $n$ sufficiently large. We have used that $\log(1+x)=x(1+o(1))$ for $x\rightarrow 0$. The last line of the equation converges to zero at a summable rate since $h_{+}\geq\log(n)^{3/d}n^{-1/(2d)}$ by assumption which concludes the proof of the uniform almost sure convergence of $\widehat{f}_{S,{\bm{\Theta}}}$. ∎

Let us turn to the standard deviation

[TABLE]

and its estimator defined in (3.13). The following lemma shows that it is uniformly bounded from above and below. The proof is deferred to Appendix B.

Lemma A.3**.**

Under Assumptions 1-3 there exist universal constants $C_{1},C_{2},n_{0}>0$ such that for any $n>n_{0},$

[TABLE]

The proof is given in Appendix B. The consistency of the estimates $\widehat{f}_{\bm{\Theta}}$ and $\widehat{f}_{S,{\bm{\Theta}}}$ shows that $\widehat{\sigma}_{\mathbf{t},h,\mathbf{v}}$ is a consistent estimator of the standard deviation $\sigma_{\mathbf{t},h,\mathbf{v}}$.

Lemma A.4**.**

Under Assumptions 1- 3,

[TABLE]

Proof.

By Lemma A.3,

[TABLE]

By Assumption 3, $f_{\bm{\Theta}}$ is uniformly bounded from below. Thus, $\widetilde{f}_{\bm{\Theta}}$ is almost surely uniformly bounded from below for sufficiently large $n$ by Lemma A.1. This shows that

[TABLE]

Here we used the boundedness of $f_{S,{\bm{\Theta}}}$ and

[TABLE]

by Lemma A.2. The claim follows now from the integrability of \big{(}(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})(s)\big{)}^{2} proved in Lemma 3.1 (ii). ∎

Lemma A.5**.**

Under Assumptions 2 resp. 2’ and 3 we have

[TABLE]

Proof.

This is a direct consequence of

[TABLE]

as shown in the proof of Lemma A.4 and

[TABLE]

∎

We discussed in Section 3 that the test statistic $T_{\mathbf{t},h,\mathbf{v}}$ relies on the unknown density $f_{\bm{\Theta}}$ and therefore we introduced the statistic $\widehat{T}_{\mathbf{t},h,\mathbf{v}}$, where the density $f_{\bm{\Theta}}$ is replaced by the estimate $\widetilde{f}_{\bm{\Theta}}$. An important part of the proof of Theorem 3.2 consists of showing that this replacement is asymptotically negligible. To this end, the bias of the estimate $1/\widetilde{f}_{\bm{\theta}}({\bm{\theta}})$ has to be controlled.

Lemma A.6**.**

Under Assumption 3,

[TABLE]

Proof.

Uniformly over ${\bm{\theta}}\in\mathbb{S}^{d-1},$

[TABLE]

Following the line of arguments in the proof of Lemma A.1, it is easy to see that \mathbb{P}\big{(}\exists{\bm{\theta}}^{\prime}:\widehat{f}_{\bm{\Theta}}({\bm{\theta}}^{\prime})\leq{f}_{\bm{\Theta}}({\bm{\theta}}^{\prime})/2\big{)} decays at a rate which is faster than polynomial. In particular, $\widetilde{f}_{\bm{\Theta}}=\widehat{f}_{\bm{\Theta}}$ except on a set with probability decaying faster than any polynomial, which concludes the proof using Lemma A.1 (i). ∎

Appendix B Proofs related to properties of the Radon transform

Proof of (3.9): Recall that in the random coefficients model with intercept

[TABLE]

where $\zeta_{1}$ is a Rademacher variable. Hence, we obtain

[TABLE]

∎

Proof of Lemma 3.1.

By assumption, $\phi_{\mathbf{t},h}$ is radially symmetric and satisfies (3.2). We fix a direction $\mathbf{v}\in\mathbb{S}^{d-1}$ and consider the directional derivative

[TABLE]

where $\phi^{\prime}$ is the usual derivative of $\phi$. The Radon transform of this directional derivative is

[TABLE]

Set $\widetilde{s}=h^{-1}(s-\langle\mathbf{t},{{\bm{\theta}}}\rangle).$ For $d>2,$ using the definition of $\widetilde{\phi}$ in (3.6),

[TABLE]

For $d=2$ let $\mathbf{w}\perp{\bm{\theta}}$ with $\|\mathbf{w}\|=1$ and write $\mathbf{b}={\bm{\theta}}\widetilde{s}+r\mathbf{w}$ for $r\in\mathbb{R}$. Then

[TABLE]

as $r\mapsto\frac{\phi^{\prime}\left(\sqrt{\widetilde{s}^{2}+r^{2}}\right)}{\sqrt{\widetilde{s}^{2}+r^{2}}}$ is an even function. Now we can proceed similarly as in the case $d>2$.

If $d$ is odd, the proof of the representation of $A(\partial_{\mathbf{v}}\phi_{\mathbf{t},h})(s,{{\bm{\theta}}})$ is completed by taking the $(d-1)$-th derivative with respect to the variable $s$. If $d$ is even, for any function $f:\mathbb{R}\rightarrow\mathbb{R}$, for any fixed $z\in\mathbb{R}$, and for $h>0$

[TABLE]

by substitution. That $\mathcal{H}_{d}f$ exists is shown below for the choice $f=\widetilde{\phi}^{(d-1)}$. Hence, we obtain $A(\partial_{\mathbf{v}}\phi_{\mathbf{t},h})(s,{\bm{\theta}})=\langle{\bm{\theta}},\mathbf{v}\rangle h^{-d-1}(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})(\tfrac{s-\langle\mathbf{t},{\bm{\theta}}\rangle}{h})$ for $d$ even.

Next we prove $\|\widetilde{\phi}^{(k)}\|_{\infty}<\infty$ for $k=0,\ldots,d+1.$ The case $k=0$ is obvious. We use the chain rule for higher order derivatives given by Faà di Bruno’s formula

[TABLE]

where $\mathcal{M}_{k}$ is the set of all $k$-tuples of non-negative integers satisfying $\sum_{j=1}^{k}jm_{j}=k$. Since $z\mapsto\phi(\sqrt{z^{2}+r^{2}})$ is a.e. $(k+1)$-times continuously differentiable, we can interchange the integral with the $k$-fold differentiation for the variable $z$ provided that

[TABLE]

exists for all $k=1,\ldots,d+1$. Applying (4.1) with $f_{1}=\phi$ and $f_{2}=\sqrt{\cdot^{2}+r^{2}}$ gives

[TABLE]

for suitable constants $C_{m_{1},\ldots,m_{k+1}}$ and $M=\sum_{j=1}^{k+1}m_{j}.$ Applying the chain rule to $f_{2}^{(j)}$ yields

[TABLE]

for non-negative integers $\ell_{j},k_{l}$ and suitable constants $C_{\ell_{j},k_{j}}$. As $\phi$ is compactly supported, it remains to show that each of the functions

[TABLE]

for $|z|\leq 1$ is uniformly bounded, where $\ell_{j},k_{j}$ are arbitrary elements of the set $\{\ell_{j},k_{j}:\ell_{j}+2k_{j}=j\},\;j=1,\ldots,k+1$. Notice that

[TABLE]

A uniform bound for the integral on the right hand side of (4.2) can be found easily when $z$ is bounded away from zero. We can thus assume that $|z|\leq\sqrt{1-z^{2}}$. Splitting the integral $\int_{0}^{\sqrt{1-z^{2}}}=\int_{0}^{|z|}+\int_{|z|}^{\sqrt{1-z^{2}}}$ and using that by Taylor expansion and Assumption 1,

[TABLE]

as well as $\max\{z^{2},r^{2}\}\leq z^{2}+r^{2}\leq 2\max\{z^{2},r^{2}\}$, we obtain an upper bound (up to some constant) for the integral on the right hand side of (4.2) by

[TABLE]

By the use of $\ell_{j}+2k_{j}=j$, $\sum_{j=1}^{k+1}jm_{j}={k+1}$ and $k\leq d+1$, we find that this is bounded by $z^{-3+M+\max\{3-M,0\}}+1$ which proves the result.

Next, we prove that $\mathcal{H}_{d}\widetilde{\phi}^{(d-1)}$ exists. Recall that $\|\widetilde{\phi}^{(d)}\|_{\infty}<\infty$ and consequently, $\widetilde{\phi}^{(d-1)}$ is Lipschitz continuous. For any Lipschitz continuous function $f$ with compact support,

[TABLE]

where $\lambda(\,\textnormal{supp}f)$ denotes the Lebesgue measure of the support of $f$. Moreover,

[TABLE]

By the Lipschitz-continuity of $f$, $|f(u-x)-f(u+x)|\lesssim|x|$ such that the r.h.s. can be bounded by a constant that does not depend on $u$. The result follows with $f=\widetilde{\phi}^{(d-1)}$. This proves assertion (i) in the Lemma.

Finally, we prove (ii). As shown above, $\widetilde{\phi}^{(d-1)}$ is bounded. For odd dimension $d$ the claim therefore follows from substitution and the compact support of $\widetilde{\phi}^{(d-1)}$. For $d$ even, substitution and the fact that the Hilbert transform $\mathcal{H}_{d}$ defines a bounded operator $L^{k}(\mathbb{R})\rightarrow L^{k}(\mathbb{R})$ for all $1<k<\infty$ yield the required result.

∎

Proof of Lemma A.3.

The existence of a uniform upper bound of $\sigma_{\mathbf{t},h,\mathbf{v}}$ follows directly from the boundedness of $f_{S,{\bm{\Theta}}}$. The uniform lower bound of $f_{\bm{\Theta}}$ follows from Assumption 3. The integrability of \big{(}\mathcal{H}_{d}(\widetilde{\phi}^{(d-1)})(s)\big{)}^{2} is shown in the proof of Lemma 3.1 (ii). For the lower bound of $\sigma_{\mathbf{t},h,\mathbf{v}}$ recall that

[TABLE]

By Assumption 2, $f_{\bm{\beta}}(\mathbf{b})\geq c_{\bm{\beta}}>0$ for all $\mathbf{b}\in[\mathbf{a}_{1},\mathbf{a}_{2}]$ and $f_{\bm{\beta}}$ is uniformly continuous. Hence, there exists $\delta>0$, which does not depend on $h$, such that $f_{\bm{\beta}}$ is uniformly bounded from below in the ball $B_{\delta}(\mathbf{t})$ of radius $\delta$ around any $\mathbf{t}\in[\mathbf{a}_{1},\mathbf{a}_{2}]$, say, $f_{{\bm{\beta}}}(\mathbf{b})>c_{{\bm{\beta}}}/2$ for all $\mathbf{b}\in\bigcup_{\mathbf{t}\in[\mathbf{a}_{1},\mathbf{a}_{2}]}B_{\delta}(\mathbf{t})$. Define for $s^{2}<\delta^{2}/(dh^{2})$

[TABLE]

where ${\bm{\theta}}_{2}^{\bot},\ldots,{\bm{\theta}}_{d}^{\bot}$ form an orthonormal basis of the orthogonal complement of $\rm span\{{\bm{\theta}}\}$. Clearly, $\mu_{d-1}(A_{\delta,\mathbf{t},h})=(2\delta)^{d-1}d^{(1-d)/2}>0$, and all $\mathbf{b}\in A_{\delta,\mathbf{t},h}$ satisfy

[TABLE]

In particular, $A_{\delta,\mathbf{t},h}\subset B_{\delta}(\mathbf{t}).$ Thus,

[TABLE]

Hence, $Rf_{\bm{\beta}}(\langle\mathbf{t},{\bm{\theta}}\rangle+hs,{\bm{\theta}})$ is uniformly bounded from below for all $\mathbf{a}_{1}+h\leq\mathbf{t}\leq\mathbf{a}_{2}-h$, ${\bm{\theta}}\in\mathbb{S}^{d-1}$ and $|s|<\delta/(\sqrt{d}h)$. Therefore,

[TABLE]

where the inequality holds uniformly over $\mathcal{T}$.

In quantum homodyne tomography, Assumption 2’ (iii) yields

[TABLE]

for $s^{2}<\delta^{2}/(dh^{2})$ if $\delta$ is sufficiently small. Hence, (4.3) holds in this case as well. Furthermore, since $\mathcal{H}_{d}(\widetilde{\phi}^{(d-1)})\in L^{2}(\mathbb{R}),$ we obtain

[TABLE]

If $\|\mathcal{H}_{d}(\widetilde{\phi}^{(d-1)})\|_{2}\neq 0$ there exists $n_{0}=n_{0}(\delta,d,\phi)\in\mathbb{N}$ such that

[TABLE]

for all $n>n_{0}$. The equality on the r.h.s. is trivial for odd dimensions $d$ and follows for even dimensions from the anti self-adjointness of the Hilbert transform and $\mathcal{H}_{d}\mathcal{H}_{d}f=-f$. ∎

Appendix C Proof of Theorem 3.2

If $\|\widetilde{\phi}^{(d-1)}\|_{2}=0$, Theorem 3.2 obviously holds. In the following we assume $\|\widetilde{\phi}^{(d-1)}\|_{2}\neq 0$ and define

[TABLE]

where the equality follows from Lemma 3.1 (i).

C.1 Controlling the effect of density estimation in the test statistic

Theorem C.1**.**

Under the assumptions of Theorem 3.2,

[TABLE]

Proof.

By the triangle inequality

[TABLE]

with $U_{\mathbf{t},h,\mathbf{v}}:=|\widehat{T}_{\mathbf{t},h,\mathbf{v}}-T_{\mathbf{t},h,\mathbf{v}}-\mathbb{E}[\widehat{T}_{\mathbf{t},h,\mathbf{v}}-T_{\mathbf{t},h,\mathbf{v}}]|$ and $V_{\mathbf{t},h,\mathbf{v}}:=|\mathbb{E}[\widehat{T}_{\mathbf{t},h,\mathbf{v}}-T_{\mathbf{t},h,\mathbf{v}}]|.$ We first bound $V_{\mathbf{t},h,\mathbf{v}}$ using

[TABLE]

and

[TABLE]

The last inequality follows for odd dimension $d$ by the boundedness of $f_{S,{\bm{\Theta}}}$ and the integrability of $\widetilde{\phi}^{(d-1)}.$ For even dimension, recall that $\mathcal{H}_{d}\widetilde{\phi}^{(d-1)}$ is bounded as shown in the proof of Lemma 3.1. Notice that

[TABLE]

by $|(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})(s)|\lesssim(1+s^{2})^{-1/2}$ (which holds for any function with compact support and bounded Hilbert transform) and the integrability of $(1+s^{2})^{-1/2}\log(s)^{-2}$ for $s\geq 2$. For the remainder, we find

[TABLE]

by the boundedness of $s\mapsto\log(|s|)^{2}f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})$ for all $|s|\geq 2,{\bm{\theta}}\in\mathbb{S}^{d-1},$ and

[TABLE]

as $\log(h/2)\geq-\log(s)/2$ for $s\geq 4/h^{2}.$ A similar argument can be used to bound the integral $\int_{-\infty}^{-2}|a_{\mathbf{t},h,\mathbf{v}}(s,{\bm{\theta}})f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})|ds.$ Applying Lemma A.6 with bandwidth $h_{*}=\log(n)^{7/(d-1)}n^{-1/(d-1)}$ gives

[TABLE]

Next, we prove $\rho_{n}:=\mathbb{P}(\sup_{(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}}U_{\mathbf{t},h,\mathbf{v}}\geq\delta_{n})\rightarrow 0$ as $n\rightarrow\infty,$ where $\delta_{n}:={(n\log(n))}^{-1/2}$. If for some positive constant $c$

[TABLE]

then by Lemma A.1

[TABLE]

for sufficiently large $c.$ Now we apply Bernstein’s inequality to

[TABLE]

By Lemma 3.1 (i), $|a_{\mathbf{t},h,\mathbf{v}}(s,{\bm{\theta}})|$ can be bounded by a constant uniformly over $(s,{\bm{\theta}})\in\mathbb{R}\times\mathbb{S}^{d-1}.$ Moreover,

[TABLE]

The inequality $\widetilde{f}_{\bm{\Theta}}\geq\log(n)^{-1}$, the uniform lower bound of $f_{\bm{\Theta}}$, $\widetilde{f}_{\bm{\Theta}}=\widehat{f}_{\bm{\Theta}}$ almost surely for $n$ sufficiently large, Lemma A.1, and the definition of $h_{*}$ imply that each summand in $U_{\mathbf{t},h,\mathbf{v}}$ is bounded on $A_{n}$ by

[TABLE]

for some constants $C,C_{1}>0$. By a change of variables in the integral for the variable $s$, the uniform boundedness of $f_{S,{\bm{\Theta}}}$, and the integrability of $a_{\mathbf{t},h,\mathbf{v}}^{2}$ as shown in Lemma 3.1 (i), we find for the conditional variance with a similar argument as above

[TABLE]

with some constants $C,C_{2}>0$. Bernstein’s inequality yields

[TABLE]

as $h_{\,\textnormal{min}}\geq n^{-1}$. Finally, the claim follows from $\beta_{h}\lesssim\frac{\sqrt{\log(n)}}{\log\log(n)}$ and the boundedness from below of $\sigma_{\mathbf{t},h,\mathbf{v}}$ shown in Lemma A.3. ∎

C.2 Approximation of the limit statistic

Define the process

[TABLE]

Note that $X_{\mathbf{t},h,\mathbf{v}}$ corresponds to the process $\widehat{X}_{\mathbf{t},h,\mathbf{v}}$ where the density estimators have been replaced by the true densities. The proof of Theorem 3.2 relies on a recently obtained Gaussian approximation result which is reproduced here for convenience.

Theorem C.2** (Chernozhukov et al., (2017), Proposition 2.1).**

Let $\mathbf{X}_{1},\ldots,\mathbf{X}_{n}$ be independent random vectors in $\mathbb{R}^{2p}$ with $\mathbb{E}[\mathrm{X}_{i,j}]=0$ and $\mathbb{E}[\mathrm{X}_{i,j}^{2}]<\infty$ for $i=1,\ldots,n,\;j=1,\ldots,2p$. Moreover, let $\mathbf{Y}_{1},\ldots,\mathbf{Y}_{n}$ be independent random vectors in $\mathbb{R}^{2p}$ with $\mathbf{Y}_{i}\sim N(\mathbf{0},\mathbb{E}[\mathbf{X}_{i}\mathbf{X}_{i}^{\top}]),\;i=1,\ldots,n$. Let $b,q>0$ be some constants and let $B_{n}\geq 1$ be a sequence of constants, possibly growing to infinity as $n\rightarrow\infty$. Denote further by $\mathcal{A}_{2p}^{\prime}$ the set of all hyperrectangles in $\mathbb{R}^{2p}$ of the form A=\big{\{}\bm{x}\in\mathbb{R}^{2p}:\bm{a}\leq\bm{x}\leq\bm{b}\big{\}} for $-\infty\leq\bm{a}\leq\bm{b}\leq\infty$. Assume that

[TABLE]

and define

[TABLE]

Then there exists a constant $C$ only depending on $b$ and $q,$ such that

[TABLE]

Theorem C.3**.**

Under the assumptions of Theorem 3.2,

[TABLE]

Proof.

To take absolute values into account, we introduce the set

[TABLE]

Moreover, for $i=1,\ldots,n,$ let $\mathbf{X}_{i}:=(\mathrm{X}_{i,1},\ldots,\mathrm{X}_{i,2p})^{\top}$ with

[TABLE]

Notice that $\sum_{i=1}^{n}X_{i,j}=n\sigma_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}^{-1}(T_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}-\mathbb{E}[T_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}])$. In a first step, we show that for $\mathbf{Z}\sim N(\mathbf{0},\mathbb{E}[\mathbf{X}_{1}\mathbf{X}_{1}^{\top}]),$

[TABLE]

Observe that by (4.1) and the uniform lower bound of $\sigma_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}$ established in Lemma A.3,

[TABLE]

Because of this bound, the expectation $\mathbb{E}[\Upsilon_{j}(S_{1},{\bm{\Theta}}_{1})]$ in the definition of $X_{i,j}$ will only provide terms of negligible order if we check the conditions of Theorem C.2. In particular, condition (i) is a direct consequence of the definition of $\sigma_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}$ in (4.6). By Lemma 3.1 (ii), the uniform lower bound of $\sigma_{\mathbf{t}_{j},h_{j},\mathbf{v}_{j}}$ in Lemma A.3, the lower bound of $f_{\bm{\Theta}}$, and the boundedness of $f_{S,{\bm{\Theta}}}$, we find for $k=1,2,$ $\max_{j=1,\ldots,2p}\mathbb{E}[|\Upsilon_{j}(S_{1},{\bm{\Theta}}_{1})|^{2+k}]\lesssim h_{\min}^{-k/2}$. This implies condition (ii) of Theorem C.2 with $B_{n}\asymp h_{\,\textnormal{min}}^{-1/2}$.

Lemma 3.1 (i) implies $\max_{j=1,\ldots,2p}|\mathrm{X}_{i,j}|\lesssim h_{\,\textnormal{min}}^{-1/2}$ which proves assertion (iii) in the theorem for any $q>0$ and $B_{n}=ch_{\,\textnormal{min}}^{-1/2}$, provided that the constant $c$ is chosen sufficiently large. Consequently, Theorem C.2 applies and for $\mathbf{Z}\sim N(\mathbf{0},\mathbb{E}[\mathbf{X}_{1}\mathbf{X}_{1}^{\top}])$

[TABLE]

choosing $q$ large enough and using Assumption 4.

In a second step, we show that there exists a version of the Gaussian noise $W$ such that

[TABLE]

To this end, we define the Gaussian process $(\widetilde{W}(f))_{f\in L^{\infty}(\mathcal{Z})}$ indexed by $L^{\infty}(\mathcal{Z})$ as the centered Gaussian process with covariance function

[TABLE]

Thus, there exists a version of $\widetilde{W}(f)$ such that \mathbf{Z}=\big{(}\widetilde{W}(\Upsilon_{1}),\ldots,\widetilde{W}(\Upsilon_{2p})\big{)}^{\top}. Recall that $(W(f))_{f\in L^{2}(\nu)}$ defines a Gaussian process whose mean and covariance functions are [math] and $\int_{\mathbb{S}^{d-1}}\int_{\mathbb{R}}f_{1}(s,{\bm{\theta}})f_{2}(s,{\bm{\theta}})dsd{\bm{\theta}}$, respectively. Basic calculations show that there exists a version of $W$ such that

[TABLE]

Hence,

[TABLE]

By (4.5), $|\int_{\mathbb{S}^{d-1}}\int_{\mathbb{R}}\Upsilon_{j}(s,{\bm{\theta}}){f_{S,{\bm{\Theta}}}(s,{\bm{\theta}})}dsd{\bm{\theta}}|\lesssim\log(h_{j})^{2}\sqrt{h_{j}}.$ Furthermore, $W(\sqrt{f_{S,{\bm{\Theta}}}})\sim N(0,1)$ which implies that $\mathbb{E}[\max_{j=1,\ldots,2p}|\widetilde{W}(\Upsilon_{j})-W(\Upsilon_{j}\sqrt{f_{S,{\bm{\Theta}}}})|]\lesssim\log(h_{\,\textnormal{max}})^{2}\sqrt{h_{\,\textnormal{max}}}.$ An application of Markov’s inequality finally proves

[TABLE]

The insertion of the bandwidth normalization terms has no influence on the convergence as translation and multiplication preserve the interval structure. ∎

C.3 Boundedness of the limit statistic

Recall from Lemma A.3 that $\sigma_{\mathbf{t},h,\mathbf{v}}$ is uniformly bounded from below whenever $h$ is sufficiently small, where the upper bound $\overline{h}$ for $h$ only depends on $\phi$, $d$ and $f_{\bm{\beta}}.$ We therefore introduce the set

[TABLE]

Theorem C.4**.**

Under the assumptions of Theorem 3.2, $\sup_{(\mathbf{t},h,\mathbf{v})\in\overline{\mathcal{T}}}\beta_{h}(|X_{\mathbf{t},h,\mathbf{v}}|/\sigma_{\mathbf{t},h,\mathbf{v}}-\alpha_{h})$ is almost surely bounded.

Proof.

We apply Theorem 6.1 in Dümbgen and Spokoiny, (2001) to the non-normalized process

[TABLE]

Denote by $\rho$ the canonical pseudo-metric on $\mathcal{T}$, induced by $Y_{\mathbf{t},h,\mathbf{v}}$

[TABLE]

In the next step, we prove

[TABLE]

By the uniform lower and upper bound for $\sigma_{\mathbf{t},h,\mathbf{v}},$

[TABLE]

In order to bound the expectation of the first term on the right hand side of (4.7), we use the boundedness properties of $f_{\bm{\Theta}}$ and $f_{S,{\bm{\Theta}}}$

[TABLE]

We show that the three terms can be bounded by the squared r.h.s. in (4.6). From Lemma 3.1 (ii) we obtain $\rho_{1}\lesssim h\|\mathbf{v}-\mathbf{v}^{\prime}\|^{2}.$ For $\rho_{2},$ we distinguish between the cases $\|\mathbf{t}-\mathbf{t}^{\prime}\|>h$ and $\|\mathbf{t}-\mathbf{t}^{\prime}\|\leq h$. In the first case, the triangle inequality and Lemma 3.1 (ii) give $\rho_{2}\lesssim h<\|\mathbf{t}-\mathbf{t}^{\prime}\|.$ In the second case, the integral w.r.t. the variable $s$ in $\rho_{2}$ is equal to

[TABLE]

Recall that the Hilbert transform and the differentiation operator commute. Therefore, using the differentiability of $\widetilde{\phi}^{(d-1)}$ which has been shown in the proof of Lemma 3.1, we find that

[TABLE]

for some $\xi$ between $s$ and $s+\tfrac{\langle\mathbf{t}-\mathbf{t}^{\prime},{\bm{\theta}}\rangle}{h}$. Hence,

[TABLE]

Here, we used the boundedness of $\mathcal{H}_{d}\widetilde{\phi}^{(d-1)}$ shown in Lemma 3.1 (i), $|\xi|\geq|s|/2$ for all $|s|\geq 2$ in the case $\|\mathbf{t}-\mathbf{t}^{\prime}\|\leq h$, and $|(\mathcal{H}_{d}\widetilde{\phi}^{(d)})(u)|\lesssim(1+u^{2})^{-1}$. The latter is obvious for $d$ odd. For $d$ even we find that $\widetilde{\phi}$ is an odd function and therefore $\widetilde{\phi}^{(d)}$ is an odd function. Moreover, for any odd function $f$ such that $\mathcal{H}_{d}f$ exists, we have, up to some constant,

[TABLE]

Here all integrals are understood in the principal value sense. Finally, a similarly argument as in the proof of Lemma 3.1 shows that $\mathcal{H}_{d}\widetilde{\phi}^{(d)}$ exists and that $|(\mathcal{H}_{d}\widetilde{\phi}^{(d)})(u)|\lesssim(1+u^{2})^{-1}$ by the compact support of $\widetilde{\phi}$.

We finally turn to $\rho_{3}$. Without loss of generality, we may assume $h\leq h^{\prime}$. We study the cases $h\leq h^{\prime}/2$ and $h>h^{\prime}/2$, separately. In the first case, the triangle inequality and Lemma 3.1 (ii) give $\rho_{3}\lesssim h+h^{\prime}\lesssim|h^{\prime}-h|.$ If $h^{\prime}/2<h\leq h^{\prime},$ we argue as for the upper bound of $\rho_{2}$ and find

[TABLE]

for some $\xi$ between $s$ and $\tfrac{h}{h^{\prime}}s$. Recall that $|(\mathcal{H}_{d}\widetilde{\phi}^{(d)})(u)|\lesssim(1+u^{2})^{-1}$ and $|(\mathcal{H}_{d}\widetilde{\phi}^{(d-1)})(u)|\lesssim(1+u^{2})^{-1/2}$. Thus,

[TABLE]

where we used that $|\xi|\geq\frac{h}{h^{\prime}}|s|>|s|/2$. Finally, $|h^{2}/h^{\prime}-h|\leq h^{\prime}-h$ implies $\rho_{3}\lesssim|h^{\prime}-h|.$

For the second term on the right hand side of (4.7) we use that $\operatorname{Var}(Y_{\mathbf{t},h,\mathbf{v}})=h$. Using again the uniform boundedness from above and below of $\sigma_{\mathbf{t},h,\mathbf{v}}$ and the fact that $|\sqrt{x}-\sqrt{y}|^{2}\leq|x-y|$ for all $x,y\geq 0$ gives

[TABLE]

For the second term in the last line, the bounds above apply which completes the proof for (4.6).

Set $\sigma^{2}(\mathbf{t},h,\mathbf{v}):=h$ and \widetilde{\rho}\bigl{(}(\mathbf{t},h,\mathbf{v}),(\mathbf{t}^{\prime},h^{\prime},\mathbf{v}^{\prime})\bigr{)}:=\big{(}\|\mathbf{v}-\mathbf{v}^{\prime}\|^{2}+{\|\mathbf{t}-\mathbf{t}^{\prime}\|}+{|h-h^{\prime}|}\big{)}^{1/2}, such that

[TABLE]

For fixed \bigl{(}(\mathbf{t},h,\mathbf{v}),(\mathbf{t}^{\prime},h^{\prime},\mathbf{v}^{\prime})\bigr{)}\in\overline{\mathcal{T}}\times\overline{\mathcal{T}}, the random variable $Y_{\mathbf{t},h,\mathbf{v}}-Y_{\mathbf{t}^{\prime},h^{\prime},\mathbf{v}^{\prime}}$ follows a normal distribution with mean zero and variance bounded by a constant multiple of \widetilde{\rho}^{2}\bigl{(}(\mathbf{t},h,\mathbf{v}),(\mathbf{t}^{\prime},h^{\prime},\mathbf{v}^{\prime})\bigr{)}. Thus, there exists a constant $M>0$ such that for any $\eta>0$,

[TABLE]

Furthermore, \mathbb{P}\bigl{(}Y_{\mathbf{t},h,\mathbf{v}}>\sqrt{h}\eta\bigr{)}\lesssim\exp(-\eta^{2}/2), as $h^{-1/2}Y_{\mathbf{t},h,\mathbf{v}}$ corresponds to a standard normal distributed random variable. Thus, conditions (i) and (ii) of Theorem 6.1 in Dümbgen and Spokoiny, (2001) are satisfied. As in Eckle et al., 2017a one shows that condition (iii) of Theorem 6.1 in Dümbgen and Spokoiny, (2001) holds with $V=(3d-1)/2$ and that the process $Y_{\mathbf{t},h,\mathbf{v}}$ is almost surely continuous on $\overline{\mathcal{T}}$ with respect to ${\rho}$. The boundedness of \sup_{(\mathbf{t},h,\mathbf{v})\in\overline{\mathcal{T}}}\bigl{(}\beta_{h}\frac{|X_{\mathbf{t},h,\mathbf{v}}|}{\sigma_{\mathbf{t},h,\mathbf{v}}}-\alpha_{h}\beta_{h}\bigr{)} follows by an application of Theorem 6.1 and Remark 1 in Dümbgen and Spokoiny, (2001). ∎

C.4 Replacing the true densities in the limit process by estimators

Theorem C.5**.**

Under the assumptions of Theorem 3.2,

[TABLE]

Proof.

Recall the definition of the symmetrized set $\mathcal{T}_{n}^{\prime}$ in (4.3) and let

[TABLE]

Lemma A.5 and an argument as in the proof of Lemma A.4 show that \sup_{(s,{\bm{\theta}})\in\mathbb{R}\times\mathbb{S}^{d-1}}\big{|}\widehat{F}(s,{\bm{\theta}})\big{|}=O(\log(n)^{-1}) for $n\rightarrow\infty$ almost surely. Define

[TABLE]

and $\Delta_{\infty,\mathbf{t},h,\mathbf{v}}:=\log(n)^{-1}X_{\mathbf{t},h,\mathbf{v}}.$ We write $\widetilde{\mathbb{P}}$ and $\widetilde{\mathbb{E}}$ for the probability and expectation conditionally on $(S_{i},{\bm{\Theta}}_{i})$, $i=n+1,\ldots,2n$. Under $\widetilde{\mathbb{P}}$, the vectors $(\Delta_{\mathbf{t},h,\mathbf{v}})_{(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}^{\prime}}$ and $(\Delta_{\infty,\mathbf{t},h,\mathbf{v}})_{(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}^{\prime}}$ are centered and normally distributed with

[TABLE]

almost surely. Hence, an application of Theorem 2.2.5 in Adler and Taylor, (2007) gives

[TABLE]

Moreover, by the almost sure asymptotic boundedness of $\sup_{(\mathbf{t},h,\mathbf{v})\in{\mathcal{T}_{n}^{\prime}}}\beta_{h}(|X_{\mathbf{t},h,\mathbf{v}}|/\sigma_{\mathbf{t},h,\mathbf{v}}-\alpha_{h})$ proved in Theorem C.4, we have $|\widetilde{\mathbb{E}}(\sup_{(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}^{\prime}}\Delta_{\infty,\mathbf{t},h,\mathbf{v}})|=O(\log(n)^{-1/2})$ almost surely. Finally, for some constant $C>0$,

[TABLE]

for $n\rightarrow\infty$ almost surely, by Markov’s inequality. The constants introduced above do not depend on the second sample $(S_{i},{\bm{\Theta}}_{i})$, $i=n+1,\ldots,2n$ and therefore the claim follows by an application of the law of iterated expectations. ∎

C.5 Replacement of the standard deviation by an estimator

Theorem C.6**.**

Under the assumptions of Theorem 3.2,

[TABLE]

Proof.

We only prove (i) as (ii) follows by a similar argument. By Lemma A.3 and Lemma A.4, $\sigma_{\mathbf{t},h,\mathbf{v}}$ and $\widehat{\sigma}_{\mathbf{t},h,\mathbf{v}}$ are almost surely uniformly bounded from below for all sufficiently large $n$. Thus,

[TABLE]

almost surely. The claim follows from Lemma A.4, Theorems C.1 and C.3 and the almost sure boundedness of $\sup_{(\mathbf{t},h,\mathbf{v})\in{\mathcal{T}_{n}}}\beta_{h}(|X_{\mathbf{t},h,\mathbf{v}}|/\sigma_{\mathbf{t},h,\mathbf{v}}-\alpha_{h})$ established in Theorem C.4. ∎

Appendix D Proofs of Theorems 3.3 and 3.4

Proof of Theorem 3.3.

We have

[TABLE]

for $n\rightarrow\infty$. Here we used (3.14) for the first equality and Theorem 3.2 for the second. ∎

Proof of Theorem 3.4.

We assume in the following that $c_{d}>0$. The case $c_{d}<0$ can be treated similarly. The following statement can be derived similarly as in the proof of Theorem 3.3 in Eckle et al., 2017a . For a null sequence $0<(\alpha_{n})_{n\in\mathbb{N}}<1$ converging sufficiently slowly and for the set $\mathcal{T}_{n}^{\prime}\subseteq\mathcal{T}_{n}$ of all triples for which the inequality

[TABLE]

is satisfied it holds that

[TABLE]

Hence, the hypotheses (3.5) are rejected simultaneously on the set of scales $\mathcal{T}_{n}^{\prime}$ with asymptotic probability one. Moreover, for a mode $\mathbf{b}_{0}$ in $(\mathbf{a}_{1},\mathbf{a}_{2})$ of $f_{\bm{\beta}}$ and any triple $(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}^{\mathbf{b}_{0}}$, one can prove that $\partial_{\mathbf{v}}f_{\bm{\beta}}(\mathbf{b})\lesssim-h$ for all $\mathbf{b}\in\,\textnormal{supp}\phi_{\mathbf{t},h}$ by following the arguments in the proof of Theorem 10 in Eckle et al., 2017b . Consequently, $\int_{\mathbb{R}^{d}}\phi_{\mathbf{t},h}(\mathbf{b})\partial_{\mathbf{v}}f_{\bm{\beta}}(\mathbf{b})\,\textnormal{d}\mathbf{b}\lesssim-h.$

As seen in Appendices C.3-C.5

[TABLE]

for $n\rightarrow\infty$ and $\sup_{(\mathbf{t},h,\mathbf{v})\in\mathcal{T}_{n}}\beta_{h}(\frac{|X_{\mathbf{t},h,\mathbf{v}}|}{{\sigma}_{\mathbf{t},h,\mathbf{v}}}-\alpha_{h})$ is finite almost surely for $n$ sufficiently large. Moreover, $\widehat{\sigma}_{\mathbf{t},h,\mathbf{v}}$ is almost surely uniformly bounded by Lemmas A.3 and A.4 for $n$ sufficiently large, such that

[TABLE]

almost surely. In order to verify (4.1), we need to pick $h$ such that $h^{d+3/2}\gtrsim(\log(n)/n)^{1/2}.$ Thus, (4.1) holds for $h\geq C\log(n)^{\frac{1}{2d+3}}n^{-\frac{1}{2d+3}}$ with some sufficiently large constant $C>0$. ∎