Search

Detecting multiple unknown objects in noisy data is a key problem in many scientific fields, such as electron microscopy imaging. A common model for the unknown objects is the linear subspace model, which assumes that the objects can be expanded in some known basis (such as the Fourier basis). In this paper, we develop an object detection algorithm that under the linear subspace model is asymptotically guaranteed to detect all objects, while controlling the family wise error rate or the false discovery rate. Numerical simulations show that the algorithm also controls the error rate with high power in the non-asymptotic regime, even in highly challenging regimes. We apply the proposed algorithm to experimental electron microscopy data set, and show that it outperforms existing standard software.

Historically, applications of RFT in fMRI have relied on assumptions of smoothness, stationarity and Gaussianity. The first two assumptions have been addressed in Part 1 of this article series. Here we address the severe non-Gaussianity of (real) fMRI data to greatly improve the performance of voxelwise RFT in fMRI group analysis. In particular, we introduce a transformation which accelerates the convergence of the Central Limit Theorem allowing us to rely on limiting Gaussianity of the test-statistic. We shall show that, when the GKF is combined with the Gaussianization transformation, we are able to accurately estimate the EEC of the excursion set of the transformed test-statistic even when the data is non-Gaussian. This allows us to drop the key assumptions of RFT inference and enables us to provide a fast approach which correctly controls the voxelwise false positive rate in fMRI. We employ a big data \cite{Eklund2016} style validation in which we process resting state data from 7000 subjects from the UK BioBank with fake task designs. We resample from this data to create realistic noise and use this to demonstrate that the error rate is correctly controlled.

Let $(M,g)$ be a Riemannian manifold. If $\mu$ is a probability measure on $M$ given by a continuous density function, one would expect the Fr\'{e}chet means of data-samples $Q=(q_1,q_2,\dots, q_N)\in M^N$, with respect to $\mu$, to behave ``generically''; e.g. the probability that the Fr\'{e}chet mean set $\mbox{FM}(Q)$ has any elements that lie in a given, positive-codimension submanifold, should be zero for any $N\geq 1$. Even this simplest instance of genericity does not seem to have been proven in the literature, except in special cases. The main result of this paper is a general, and stronger, genericity property: given i.i.d. absolutely continuous $M$-valued random variables $X_1,\dots, X_N$, and a subset $A\subset M$ of volume-measure zero, $\mbox{Pr}\left\{\mbox{FM}(\{X_1,\dots,X_N\})\subset M\backslash A\right\}=1.$ We also establish a companion theorem for equivariant Fr\'{e}chet means, defined when $(M,g)$ arises as the quotient of a Riemannian manifold $(\widetilde{M},\tilde{g})$ by a free, isometric action of a finite group. The equivariant Fr\'{e}chet means lie in $\widetilde{M}$, but, as we show, project down to the ordinary Fr\'{e}chet sample means, and enjoy a similar genericity property. Both these theorems are proven as consequences of a purely geometric (and quite general) result that constitutes the core mathematics in this paper: If $A\subset M$ has volume zero in $M$ , then the set $\{Q\in M^N : \mbox{FM}(Q) \cap A\neq\emptyset\}$ has volume zero in $M^N$. We conclude the paper with an application to partial scaling-rotation means, a type of mean for symmetric positive-definite matrices., Comment: 36 pages, 2 figures

We study extensions of Fr\'{e}chet means for random objects in the space ${\rm Sym}^+(p)$ of $p \times p$ symmetric positive-definite matrices using the scaling-rotation geometric framework introduced by Jung et al. [\textit{SIAM J. Matrix. Anal. Appl.} \textbf{36} (2015) 1180-1201]. The scaling-rotation framework is designed to enjoy a clearer interpretation of the changes in random ellipsoids in terms of scaling and rotation. In this work, we formally define the \emph{scaling-rotation (SR) mean set} to be the set of Fr\'{e}chet means in ${\rm Sym}^+(p)$ with respect to the scaling-rotation distance. Since computing such means requires a difficult optimization, we also define the \emph{partial scaling-rotation (PSR) mean set} lying on the space of eigen-decompositions as a proxy for the SR mean set. The PSR mean set is easier to compute and its projection to ${\rm Sym}^+(p)$ often coincides with SR mean set. Minimal conditions are required to ensure that the mean sets are non-empty. Because eigen-decompositions are never unique, neither are PSR means, but we give sufficient conditions for the sample PSR mean to be unique up to the action of a certain finite group. We also establish strong consistency of the sample PSR means as estimators of the population PSR mean set, and a central limit theorem. In an application to multivariate tensor-based morphometry, we demonstrate that a two-group test using the proposed PSR means can have greater power than the two-group test using the usual affine-invariant geometric framework for symmetric positive-definite matrices.

We study asymptotic statistical inference in the space of bounded functions endowed with the supremums norm over an arbitrary metric space $S$ using a novel concept: Simultaneous COnfidence Region of Excursion (SCoRE) Sets. They simultaneously quantify the uncertainty of several lower and upper excursion sets of a target function. We investigate their connection to multiple hypothesis tests controlling the familywise error rate in the strong sense and show that they grant a unifying perspective on several statistical inference tools such as simultaneous confidence bands, quantification of uncertainties in level set estimation, for example, CoPE sets, and multiple hypothesis testing over $S$, for example, finding relevant differences or regions of equivalence within $S$. In particular, our abstract setting allows us to refine and reduce the assumptions in recent articles on CoPE sets and relevance and equivalence testing using the supremums norm.

We study power approximation formulas for peak detection using Gaussian random field theory. The approximation, based on the expected number of local maxima above the threshold $u$, $\mathbb{E}[M_u]$, is proved to work well under three asymptotic scenarios: small domain, large threshold, and sharp signal. An adjusted version of $\mathbb{E}[M_u]$ is also proposed to improve accuracy when the expected number of local maxima $\mathbb{E}[M_{-\infty}]$ exceeds 1. Cheng and Schwartzman (2018) developed explicit formulas for $\mathbb{E}[M_u]$ of smooth isotropic Gaussian random fields with zero mean. In this paper, these formulas are extended to allow for rotational symmetric mean functions, so that they are suitable for power calculations. We also apply our formulas to 2D and 3D simulated datasets, and the 3D data is induced by a group analysis of fMRI data from the Human Connectome Project to measure performance in a realistic setting.

Motivated by the questions of risk assessment in climatology (temperature change in North America) and medicine (impact of statin usage and COVID-19 on hospitalized patients), we address the problem of estimating the set in the domain of a function whose image equals a predefined subset. Existing methods that construct confidence sets require strict assumptions. We generalize the estimation of such sets to dense and non-dense domains with protection against "data peeking" by proving that confidence sets of multiple levels can be simultaneously constructed with the desired confidence non-asymptotically through inverting simultaneous confidence bands. A non-parametric bootstrap algorithm and code are provided.

The analysis of excursion sets in imaging data is essential to a wide range of scientific disciplines such as neuroimaging, climatology and cosmology. Despite growing literature, there is little published concerning the comparison of processes that have been sampled across the same spatial region but which reflect different study conditions. Given a set of asymptotically Gaussian random fields, each corresponding to a sample acquired for a different study condition, this work aims to provide confidence statements about the intersection, or union, of the excursion sets across all fields. Such spatial regions are of natural interest as they directly correspond to the questions "all random fields exceed a predetermined threshold?", or "Where does at least one random field exceed a predetermined threshold?". To assess the degree of spatial variability present, we develop a method that provides, with a desired confidence, subsets and supersets of spatial regions defined by logical conjunctions (i.e. set intersections) or disjunctions (i.e. set unions), without any assumption on the dependence between the different fields. The method is verified by extensive simulations and demonstrated using a task-fMRI dataset to identify brain regions with activation common to four variants of a working memory task., Comment: For Supplementary Theory see https://drive.google.com/file/d/1hXhrstxlHY_MMjfwE3VxPQYxp0N6JMp2/view?usp=sharing . For Supplementary Results see https://drive.google.com/file/d/156EUIYq1YIJblXI4etYmT1JRjdVpQNhf/view?usp=sharing . For code see https://github.com/TomMaullin/ConfSets

Given a functional central limit (fCLT) for an estimator and a parameter transformation, we construct random processes, called functional delta residuals, which asymptotically have the same covariance structure as the limit process of the functional delta method. An explicit construction of these residuals for transformations of moment-based estimators and a multiplier bootstrap fCLT for the resulting functional delta residuals are proven. The latter is used to consistently estimate the quantiles of the maximum of the limit process of the functional delta method in order to construct asymptotically valid simultaneous confidence bands for the transformed functional parameters. Performance of the coverage rate of the developed construction, applied to functional versions of Cohen's d, skewness and kurtosis, is illustrated in simulations and their application to test Gaussianity is discussed.

Let $\{X(t), t\in M\}$ and $\{Z(t'), t'\in M'\}$ be smooth Gaussian random fields parameterized on Riemannian manifolds $M$ and $M'$, respectively, such that $X(t) = Z(f(t))$, where $f: M \to M'$ is a diffeomorphic transformation. We study the expected number and height distribution of the critical points of $X$ in connection with those of $Z$. As an important case, when $X$ is an anisotropic Gaussian random field, then we show that its expected number of critical points becomes proportional to that of an isotropic field $Z$, while the height distribution remains the same as that of $Z$.

The expected Euler characteristic (EEC) curve of excursion sets of a Gaussian random field is used to approximate the distribution of its supremum for high thresholds. Viewed as a function of the excursion threshold, the EEC is expressed by the Gaussian kinematic formula (GKF) as a linear function of the Lipschitz-Killing curvatures (LKCs) of the field, which solely depend on the domain and covariance function of the field. So far its use for non-stationary Gaussian fields over non-trivial domains has been limited because in this case the LKCs are difficult to estimate. In this paper, consistent estimators of the LKCs are proposed as linear projections of "pinned" observed Euler characteristic curves and a linear parametric estimator of the EEC curve is obtained, which is more efficient than its nonparametric counterpart for repeated observations. A multiplier bootstrap modification reduces the variance of the estimator, and allows estimation of LKCs and EEC of the limiting field of non-Gaussian fields satisfying a functional CLT. The proposed methods are evaluated using simulations of 2D fields and illustrated in thresholding of 3D fMRI brain activation maps and cosmological simulations on the 2-sphere.

We discuss a new procedure to search for point sources in Cosmic Microwave background maps; in particular, we aim at controlling the so-called False Discovery Rate, which is defined as the expected value of false discoveries among pixels which are labelled as contaminated by point sources. We exploit a procedure called STEM, which is based on the following four steps: 1) needlet filtering of the observed CMB maps, to improve the signal to noise ratio; 2) selection of candidate peaks, i.e., the local maxima of filtered maps; 3) computation of \emph{p-}values for local maxima; 4) implementation of the multiple testing procedure, by means of the so-called Benjamini-Hochberg method. Our procedures are also implemented on the latest release of Planck CMB maps.

This article constructs simultaneous confidence bands (SCBs) for functional parameters using the Gaussian Kinematic formula of $t$-processes (tGKF). Although the tGKF relies on Gaussianity, we show that a central limit theorem (CLT) for the parameter of interest is enough to obtain asymptotically precise covering rates even for non-Gaussian processes. As a proof of concept we study the functional signal-plus-noise model and derive a CLT for an estimator of the Lipschitz-Killing curvatures, the only data dependent quantities in the tGKF SCBs. Extensions to discrete sampling with additive observation noise are discussed using scale space ideas from regression analysis. Here we provide sufficient conditions on the processes and kernels to obtain convergence of the functional scale space surface. The theoretical work is accompanied by a simulation study comparing different methods to construct SCBs for the population mean. We show that the tGKF works well even for small sample sizes and only a Rademacher multiplier-$t$ bootstrap performs similarily well. For larger sample sizes the tGKF often outperforms the bootstrap methods and is computational faster. We apply the method to diffusion tensor imaging (DTI) fibers using a scale space approach for the difference of population means. R code is available in our Rpackage SCBfda.

We study the topology generated by the temperature fluctuations of the Cosmic Microwave Background (CMB) radiation, as quantified by the number of components and holes, formally given by the Betti numbers, in the growing excursion sets. We compare CMB maps observed by the Planck satellite with a thousand simulated maps generated according to the LCDM paradigm with Gaussian distributed fluctuations. The survey of the CMB over $\mathbb{S}^2$ is incomplete due to obfuscation effects by bright point sources and other extended foreground objects like our own galaxy. To deal with such situations, where analysis in the presence of "masks" is of importance, we introduce the concept of relative homology. The parametric $\chi^2$-test shows differences between observations and simulations, yielding $p$-values at per-cent to less than per-mil levels roughly between 2 to 7 degrees. The highest observed deviation for $b_0$ and $b_1$ is approximately between $3\sigma$-4$\sigma$ at scales of 3 to 7 degrees. There are reports of mildly unusual behaviour of the Euler characteristic at 3.66 degrees in the literature, computed from independent measurements of the CMB temperature fluctuations by Planck's predecessor WMAP satellite. The mildly anomalous behaviour of Euler characteristic is related to the strongly anomalous behaviour of components and holes. These are also the scales at which the observed maps exhibit low variance compared to the simulations. Non-parametric tests show even stronger differences at almost all scales. Regardless, beyond the trivial possibility that this may still be a manifestation of an extreme Gaussian case, these observations, along with the super-horizon scales involved, may motivate to look at primordial non-Gaussianity. Alternative scenarios worth exploring may be models with non-trivial topology., Comment: 18 pages, 16 figures, 3 tables, submitted to Astronomy \& Astrophysics

In brain oncology, it is routine to evaluate the progress or remission of the disease based on the differences between a pre-treatment and a post-treatment Positron Emission Tomography (PET) scan. Background adjustment is necessary to reduce confounding by tissue-dependent changes not related to the disease. When modeling the voxel intensities for the two scans as a bivariate Gaussian mixture, background adjustment translates into standardizing the mixture at each voxel, while tumor lesions present themselves as outliers to be detected. In this paper, we address the question of how to standardize the mixture to a standard multivariate normal distribution, so that the outliers (i.e., tumor lesions) can be detected using a statistical test. We show theoretically and numerically that the tail distribution of the standardized scores is favorably close to standard normal in a wide range of scenarios while being conservative at the tails, validating voxelwise hypothesis testing based on standardized scores. To address standardization in spatially heterogeneous image data, we propose a spatial and robust multivariate expectation-maximization (EM) algorithm, where prior class membership probabilities are provided by transformation of spatial probability template maps and the estimation of the class mean and covariances are robust to outliers. Simulations in both univariate and bivariate cases suggest that standardized scores with soft assignment have tail probabilities that are either very close to or more conservative than standard normal. The proposed methods are applied to a real data set from a PET phantom experiment, yet they are generic and can be used in other contexts.

Jung et al. (2015) introduced a geometric structure on ${\rm Sym}^+(p)$, the set of $p \times p$ symmetric positive-definite matrices, based on eigen-decomposition. Eigenstructure determines both a stratification of ${\rm Sym}^+(p)$, defined by eigenvalue multiplicities, and fibers of the "eigen-composition" map $F:M(p):=SO(p)\times{\rm Diag}^+(p)\to{\rm Sym}^+(p)$. When $M(p)$ is equipped with a suitable Riemannian metric, the fiber structure leads to notions of scaling-rotation distance between $X,Y\in {\rm Sym}^+(p)$, the distance in $M(p)$ between fibers $F^{-1}(X)$ and $F^{-1}(Y)$, and minimal smooth scaling-rotation (MSSR) curves, images in ${\rm Sym}^+(p)$ of minimal-length geodesics connecting two fibers. In this paper we study the geometry of the triple $(M(p),F,{\rm Sym}^+(p))$, focusing on some basic questions: For which $X,Y$ is there a unique MSSR curve from $X$ to $Y$? More generally, what is the set ${\cal M}(X,Y)$ of MSSR curves from $X$ to $Y$? This set is influenced by two potential types of non-uniqueness. We translate the question of whether the second type can occur into a question about the geometry of Grassmannians $G_m({\bf R}^p)$, with $m$ even, that we answer for $p\leq 4$ and $p\geq 11$. Our method of proof also yields an interesting half-angle formula concerning principal angles between subspaces of ${\bf R}^p$ whose dimensions may or may not be equal. The general-$p$ results concerning MSSR curves and scaling-rotation distance that we establish here underpin the explicit $p=3$ results in Groisser et al. (2017). Addressing the uniqueness-related questions requires a thorough understanding of the fiber structure of $M(p)$, which we also provide., Comment: arXiv admin note: text overlap with arXiv:1602.01187

We present a topological multiple testing scheme for detecting peaks on the sphere under isotropic Gaussian noise, where tests are performed at local maxima of the observed field filtered by the spherical needlet transform. Our setting is different from the standard Euclidean/large same asymptotic framework, yet highly relevant to realistic experimental circumstances for some important areas of application in astronomy. More precisely, we focus on cases where a single realization of a smooth isotropic Gaussian random field on the sphere is observed, and a number of well-localized signals are superimposed on such background field. The proposed algorithms, combined with the Benjamini-Hochberg procedure for thresholding p-values, provide asymptotic strong control of the False Discovery Rate (FDR) and power consistency as the signal strength and the frequency of the needlet transform get large. This novel multiple testing method is illustrated in a simulation of point-source detection in Cosmic Microwave Background radiation (CMB) data.

We investigate a geometric computational framework, called the "scaling-rotation framework", on ${\rm Sym}^+(p)$, the set of $p \times p$ symmetric positive-definite (SPD) matrices. The purpose of our study is to lay geometric foundations for statistical analysis of SPD matrices, in situations in which eigenstructure is of fundamental importance, for example diffusion-tensor imaging (DTI). Eigen-decomposition, upon which the scaling-rotation framework is based, determines both a stratification of ${\rm Sym}^+(p)$, defined by eigenvalue multiplicities, and fibers of the "eigen-composition" map $SO(p)\times{\rm Diag}^+(p)\to{\rm Sym}^+(p)$. This leads to the notion of scaling-rotation distance [Jung et al. (2015)], a measure of the minimal amount of scaling and rotation needed to transform an SPD matrix, $X,$ into another, $Y,$ by a smooth curve in ${\rm Sym}^+(p)$. Our main goal in this paper is the systematic characterization and analysis of minimal smooth scaling-rotation (MSSR) curves, images in ${\rm Sym}^+(p)$ of minimal-length geodesics connecting two fibers in the "upstairs" space $SO(p)\times{\rm Diag}^+(p)$. The length of such a geodesic connecting the fibers over $X$ and $Y$ is what we define to be the scaling-rotation distance from $X$ to $Y.$ For the important low-dimensional case $p = 3$ (the home of DTI), we find new explicit formulas for MSSR curves and for the scaling-rotation distance, and identify ${\cal M}(X,Y)$ in all "nontrivial" cases. The quaternionic representation of $SO(3)$ is used in these computations. We also provide closed-form expressions for scaling-rotation distance and MSSR curves for the case $p = 2$., Comment: 15 figures

We obtain formulae for the expected number and height distribution of critical points of smooth isotropic Gaussian random fields parameterized on Euclidean space or spheres of arbitrary dimension. The results hold in general in the sense that there are no restrictions on the covariance function of the field except for smoothness and isotropy. The results are based on a characterization of the distribution of the Hessian of the Gaussian field by means of the family of Gaussian orthogonally invariant (GOI) matrices, of which the Gaussian orthogonal ensemble (GOE) is a special case. The obtained formulae depend on the covariance function only through a single parameter (Euclidean space) or two parameters (sphere), and include the special boundary case of random Laplacian eigenfunctions., Comment: 1 figure

A new approach to detect change points based on differential smoothing and multiple testing is presented for long data sequences modeled as piecewise constant functions plus stationary ergodic Gaussian noise. As an application of the STEM algorithm for peak detection developed in \citet{schwartzman2011multiple} and \citet{cheng2017multiple}, the method detects change points as significant local maxima and minima after smoothing and differentiating the observed sequence. The algorithm, combined with the Benjamini-Hochberg procedure for thresholding p-values, provides asymptotic strong control of the False Discovery Rate (FDR) and power consistency, as the length of the sequence and the size of the jumps get large. Simulations show that FDR levels are maintained in non-asymptotic conditions and guide the choice of smoothing bandwidth. The methods are illustrated in magnetometer sensor data and genomic array-CGH data. An R package named "dSTEM" is available in R cran.

The explicit formulae for the height distribution and expected number of local maxima have been obtained for isotropic Gaussian random fields on certain low-dimensional Euclidean space or low-dimensional spheres., Comment: 2 figures

The goal of this paper is to give confidence regions for the excursion set of a spatial function above a given threshold from repeated noisy observations on a fine grid of fixed locations. Given an asymptotically Gaussian estimator of the target function, a pair of data-dependent nested excursion sets are constructed that are sub- and super-sets of the true excursion set, respectively, with a desired confidence. Asymptotic coverage probabilities are determined via a multiplier bootstrap method, not requiring Gaussianity of the original data nor stationarity or smoothness of the limiting Gaussian field. The method is used to determine regions in North America where the mean summer and winter temperatures are expected to increase by mid 21st century by more than 2 degrees Celsius.

This article gives a formal definition of a lognormal family of probability distributions on the set of symmetric positive definite (PD) matrices, seen as a matrix-variate extension of the univariate lognormal family of distributions. Two forms of this distribution are obtained as the large sample limiting distribution via the central limit theorem of two types of geometric averages of i.i.d. PD matrices: the log-Euclidean average and the canonical geometric average. These averages correspond to two different geometries imposed on the set of PD matrices. The limiting distributions of these averages are used to provide large-sample confidence regions for the corresponding population means. The methods are illustrated on a voxelwise analysis of diffusion tensor imaging data, permitting a comparison between the various average types from the point of view of their sampling variability., Comment: 28 pages, 8 figures

We introduce a new geometric framework for the set of symmetric positive-definite (SPD) matrices, aimed to characterize deformations of SPD matrices by individual scaling of eigenvalues and rotation of eigenvectors of the SPD matrices. To characterize the deformation, the eigenvalue-eigenvector decomposition is used to find alternative representations of SPD matrices, and to form a Riemannian manifold so that scaling and rotations of SPD matrices are captured by geodesics on this manifold. The problems of non-unique eigen-decompositions and eigenvalue multiplicities are addressed by finding minimal-length geodesics, which gives rise to a distance and an interpolation method for SPD matrices. Computational procedures to evaluate the minimal scaling--rotation deformations and distances are provided for the most useful cases of $2 \times 2$ and $3 \times 3$ SPD matrices. In the new geometric framework, minimal scaling--rotation curves interpolate eigenvalues at constant logarithmic rate, and eigenvectors at constant angular rate. In the context of diffusion tensor imaging, this results in better behavior of the trace, determinant and fractional anisotropy of interpolated SPD matrices in typical cases.

A topological multiple testing scheme is presented for detecting peaks in images under stationary ergodic Gaussian noise, where tests are performed at local maxima of the smoothed observed signals. The procedure generalizes the one-dimensional scheme of Schwartzman et al. (2011) to Euclidean domains of arbitrary dimension. Two methods are developed according to two different ways of computing p-values: (i) using the exact distribution of the height of local maxima (Cheng and Schwartzman, 2014), available explicitly when the noise field is isotropic; (ii) using an approximation to the overshoot distribution of local maxima above a pre-threshold (Cheng and Schwartzman, 2014), applicable when the exact distribution is unknown, such as when the stationary noise field is non-isotropic. The algorithms, combined with the Benjamini-Hochberg procedure for thresholding p-values, provide asymptotic strong control of the False Discovery Rate (FDR) and power consistency, with specific rates, as the search space and signal strength get large. The optimal smoothing bandwidth and optimal pre-threshold are obtained to achieve maximum power. Simulations show that FDR levels are maintained in non-asymptotic conditions. The methods are illustrated in a nanoscopy image analysis problem of detecting fluorescent molecules against the image background., Comment: 30 pages, 5 figures. arXiv admin note: text overlap with arXiv:1203.3063

Let $\{f(t): t\in T\}$ be a smooth Gaussian random field over a parameter space $T$, where $T$ may be a subset of Euclidean space or, more generally, a Riemannian manifold. For any local maximum of $f(t)$ located at $t_0$ in the interior of $T$, we provide general formulae and asymptotic approximations for both the tail distribution of the height of a local maximum $\mathbb{P}\{f(t_0)>u | t_0 \text{is a local maximum of} f(t) \}$ and the overshoot distribution of a local maximum $\mathbb{P}\{f(t_0)>u+v | t_0 \text{is a local maximum of} f(t) \text{and} f(t_0)>v\}$. Assuming further that $f$ is isotropic, we apply techniques from random matrix theory related to the Gaussian orthogonal ensemble to compute such conditional probabilities explicitly when $T$ is Euclidean or a sphere of arbitrary dimension. Such calculations are motivated by the statistical problem of detecting peaks in the presence of smooth Gaussian noise., Comment: 40 pages, 2 figures

A topological multiple testing approach to peak detection is proposed for the problem of detecting transcription factor binding sites in ChIP-Seq data. After kernel smoothing of the tag counts over the genome, the presence of a peak is tested at each observed local maximum, followed by multiple testing correction at the desired false discovery rate level. Valid p-values for candidate peaks are computed via Monte Carlo simulations of smoothed Poisson sequences, whose background Poisson rates are obtained via linear regression from a Control sample at two different scales. The proposed method identifies nearby binding sites that other methods do not., Comment: Published in at http://dx.doi.org/10.1214/12-AOAS594 the Annals of Applied Statistics (http://www.imstat.org/aoas/) by the Institute of Mathematical Statistics (http://www.imstat.org)

A topological multiple testing scheme for one-dimensional domains is proposed where, rather than testing every spatial or temporal location for the presence of a signal, tests are performed only at the local maxima of the smoothed observed sequence. Assuming unimodal true peaks with finite support and Gaussian stationary ergodic noise, it is shown that the algorithm with Bonferroni or Benjamini--Hochberg correction provides asymptotic strong control of the family wise error rate and false discovery rate, and is power consistent, as the search space and the signal strength get large, where the search space may grow exponentially faster than the signal strength. Simulations show that error levels are maintained for nonasymptotic conditions, and that power is maximized when the smoothing kernel is close in shape and bandwidth to the signal peaks, akin to the matched filter theorem in signal processing. The methods are illustrated in an analysis of electrical recordings of neuronal cell activity., Comment: Published in at http://dx.doi.org/10.1214/11-AOS943 the Annals of Statistics (http://www.imstat.org/aos/) by the Institute of Mathematical Statistics (http://www.imstat.org)

This paper considers the problem of detecting equal-shaped non-overlapping unimodal peaks in the presence of Gaussian ergodic stationary noise, where the number, location and heights of the peaks are unknown. A multiple testing approach is proposed in which, after kernel smoothing, the presence of a peak is tested at each observed local maximum. The procedure provides strong control of the family wise error rate and the false discovery rate asymptotically as both the signal-to-noise ratio (SNR) and the search space get large, where the search space may grow exponentially as a function of SNR. Simulations assuming a Gaussian peak shape and a Gaussian autocorrelation function show that desired error levels are achieved for relatively low SNR and are robust to partial peak overlap. Simulations also show that detection power is maximized when the smoothing bandwidth is close to the bandwidth of the signal peaks, akin to the well-known matched filter theorem in signal processing. The procedure is illustrated in an analysis of electrical recordings of neuronal cell activity., Comment: 37 pages, 8 figures

In large scale multiple testing, the use of an empirical null distribution rather than the theoretical null distribution can be critical for correct inference. This paper proposes a ``mode matching'' method for fitting an empirical null when the theoretical null belongs to any exponential family. Based on the central matching method for $z$-scores, mode matching estimates the null density by fitting an appropriate exponential family to the histogram of the test statistics by Poisson regression in a region surrounding the mode. The empirical null estimate is then used to estimate local and tail false discovery rate (FDR) for inference. Delta-method covariance formulas and approximate asymptotic bias formulas are provided, as well as simulation studies of the effect of the tuning parameters of the procedure on the bias-variance trade-off. The standard FDR estimates are found to be biased down at the far tails. Correlation between test statistics is taken into account in the covariance estimates, providing a generalization of Efron's ``wing function'' for exponential families. Applications with $\chi^2$ statistics are shown in a family-based genome-wide association study from the Framingham Heart Study and an anatomical brain imaging study of dyslexia in children., Comment: Published in at http://dx.doi.org/10.1214/08-AOAS184 the Annals of Applied Statistics (http://www.imstat.org/aoas/) by the Institute of Mathematical Statistics (http://www.imstat.org)

This article presents maximum likelihood estimators (MLEs) and log-likelihood ratio (LLR) tests for the eigenvalues and eigenvectors of Gaussian random symmetric matrices of arbitrary dimension, where the observations are independent repeated samples from one or two populations. These inference problems are relevant in the analysis of diffusion tensor imaging data and polarized cosmic background radiation data, where the observations are, respectively, $3\times3$ and $2\times2$ symmetric positive definite matrices. The parameter sets involved in the inference problems for eigenvalues and eigenvectors are subsets of Euclidean space that are either affine subspaces, embedded submanifolds that are invariant under orthogonal transformations or polyhedral convex cones. We show that for a class of sets that includes the ones considered in this paper, the MLEs of the mean parameter do not depend on the covariance parameters if and only if the covariance structure is orthogonally invariant. Closed-form expressions for the MLEs and the associated LLRs are derived for this covariance structure., Comment: Published in at http://dx.doi.org/10.1214/08-AOS628 the Annals of Statistics (http://www.imstat.org/aos/) by the Institute of Mathematical Statistics (http://www.imstat.org)