Search

Because haplotype information is of widespread interest in biomedical applications, effort has been put into their reconstruction. Here, we propose an efficient method, called haploSep, that is able to accurately infer major haplotypes and their frequencies just from multiple samples of allele frequency data. Even the accuracy of experimentally obtained allele frequencies can be improved by re-estimating them from our reconstructed haplotypes. From a methodological point of view, we model our problem as a multivariate regression problem where both the design matrix and the coefficient matrix are unknown. Compared to other methods, haploSep is very fast, with linear computational complexity in the haplotype length. We illustrate our method on simulated and real data focusing on experimental evolution and microbial data. haploSep is a computationally efficient method to infer major haplotypes and their frequencies from multiple samples of allele frequency data, and to provide improved estimates of experimentally obtained allele frequencies.

Dilated cardiomyopathy (DCM) is a major risk factor for heart failure and is associated with the development of life-threatening cardiac arrhythmias. Using a patient-specific induced pluripotent stem cell-derived cardiomyocyte (iPSC-CM) model harbouring a mutation in cardiac troponin T (R173W), we aim to examine the cellular basis of arrhythmogenesis in DCM patients with this mutation. iPSC from control (Ctrl) and DCM-TnT-R173W donors from the same family were differentiated into iPSC-CM and analysed through optical action potential (AP) recordings, simultaneous measurement of cytosolic calcium concentration ([Ca2+]i) and membrane currents and separately assayed using field stimulation to detect the threshold for AP- and [Ca2+]i-alternans development. AP duration was unaltered in TnT-R173W iPSC-CM. Nevertheless, TnT-R173W iPSC-CM showed a strikingly low stimulation threshold for AP- and [Ca2+]i-alternans. Myofilaments are known to play a role as intracellular Ca2+ buffers and here we show increased Ca2+ affinity of intracellular buffers in TnT-R173W cells, indicating increased myofilament sensitivity to Ca2+. Similarly, EMD57033, a myofilament Ca2+ sensitiser, replicated the abnormal [Ca2+]i dynamics observed in TnT-R173W samples and lowered the threshold for alternans development. In contrast, application of a Ca2+ desensitiser (blebbistatin) to TnT-R173W iPSC-CM was able to phenotypically rescue Ca2+ dynamics, normalising Ca2+ transient profile and minimising the occurrence of Ca2+ alternans at physiological frequencies. This finding suggests that increased Ca2+ buffering likely plays a major arrhythmogenic role in patients with DCM, specifically in those with mutations in cardiac troponin T. In addition, we propose that modulation of myofilament Ca2+ sensitivity could be an effective anti-arrhythmic target for pharmacological management of this disease.

In this discussion, we compare the choice of seeded intervals and that of random intervals for change point segmentation from practical, statistical and computational perspectives. Furthermore, we investigate a novel estimator of the noise level, which improves many existing model selection procedures (including the steepest drop to low levels), particularly for challenging frequent change point scenarios with low signal-to-noise ratios., To appear in the Journal of the Korean Statistical Society

Many modern statistically efficient methods come with tremendous computational challenges, often leading to large-scale optimisation problems. In this work, we examine such computational issues for recently developed estimation methods in nonparametric regression with a specific view on image denoising. We consider in particular certain variational multiscale estimators which are statistically optimal in minimax sense, yet computationally intensive. Such an estimator is computed as the minimiser of a smoothness functional (e.g., TV norm) over the class of all estimators such that none of its coefficients with respect to a given multiscale dictionary is statistically significant. The so obtained multiscale Nemirowski-Dantzig estimator (MIND) can incorporate any convex smoothness functional and combine it with a proper dictionary including wavelets, curvelets and shearlets. The computation of MIND in general requires to solve a high-dimensional constrained convex optimisation problem with a specific structure of the constraints induced by the statistical multiscale testing criterion. To solve this explicitly, we discuss three different algorithmic approaches: the Chambolle-Pock, ADMM and semismooth Newton algorithms. Algorithmic details and an explicit implementation is presented and the solutions are then compared numerically in a simulation study and on various test images. We thereby recommend the Chambolle-Pock algorithm in most cases for its fast convergence. We stress that our analysis can also be transferred to signal recovery and other denoising problems to recover more general objects whenever it is possible to borrow statistical strength from data patches of similar object structure., Codes are available at https://github.com/housenli/MIND

Since haplotype information is of widespread interest in biomedical applications, effort has been put into their reconstruction. Here, we propose a new, computationally efficient method, called haploSep, that is able to accurately infer major haplotypes and their frequencies just from multiple samples of allele frequency data. Our approach seems to be the first that is able to estimate more than one haplotype given such data. Even the accuracy of experimentally obtained allele frequencies can be improved by re-estimating them from our reconstructed haplotypes. From a methodological point of view, we model our problem as a multivariate regression problem where both the design matrix and the coefficient matrix are unknown. The design matrix, with 0/1 entries, models haplotypes and the columns of the coefficient matrix represent the frequencies of haplotypes, which are non-negative and sum up to one. We illustrate our method on simulated and real data focusing on experimental evolution and microbial data.

Nous considerons le probleme inverse stochastique de reconstruire $f$ a partir de donnees bruitees $Y=Tf+\sigma \xi $ ou $\xi $ est un bruite blanc et $T$ un operateur compact entre espaces de Hilbert. Considerant des methodes de reconstruction generales de la forme $\hat{f}_{\alpha }=q_{\alpha }(T^{*}T)T^{*}Y$ avec un filtre ordonne $q_{\alpha }$, nous examinons le choix du parametre de regularisation $\alpha $ en minimisant un estimateur non biaise du risque predictif $\mathbb{E}[\Vert Tf-T\hat{f}_{\alpha }\Vert^{2}]$. Le parametre correspondant $\alpha_{\mathrm{pred}}$ et son utilisation sont bien connus dans la litterature mais les inegalites oracles et les resultats d’optimalite dans ce cadre general sont inconnus. Nous prouvons une inegalite oracle (generalisee), qui relie le risque direct $\mathbb{E}[\Vert f-\hat{f}_{\alpha_{\mathrm{pred}}}\Vert^{2}]$ au risque de l’oracle de prediction $\inf_{\alpha >0}\mathbb{E}[\Vert Tf-T\hat{f}_{\alpha }\Vert^{2}]$. A partir de cette inegalite oracle nous sommes alors capable de conclure que la regle de choix du parametre examine est d’ordre optimale au sens minimax. Finalement nous presentons aussi des simulations numeriques qui confirment l’optimalite de l’ordre de la methode et la qualite du choix du parametre dans des situations discretes.

In this chapter we discuss how to obtain statistical guarantees in photonic imaging. We start with an introduction to hypothesis testing in the context of imaging, more precisely we describe how to test if there is signal in a specific region of interest (RoI) or just noise. Afterwards we extend this approach to a family of RoIs and examine the occurring problems such as inflation of type I error and dependency issues. We discuss how to control the family-wise error rate by different modifications, and provide a connection to extreme value theory. Afterwards we present possible extension to inverse problems. Moving from testing to estimation, we finally introduce a method which constructs an estimator of the desired quantity of interest with automatic smoothness guarantees.

Many multiscale segmentation methods have been proven to work successfully for detecting multiple change-points, mainly because they provide faithful statistical statements, while at the same time allowing for efficient computation. Underpinning theory has been studied exclusively for models which assume that the signal is an unknown step function. However, when the signal is only approximately piecewise constant, which often occurs in practical applications, the behavior of multiscale segmentation methods is still not well studied. To narrow this gap, we investigate the asymptotic properties of a certain class of \emph{multiscale change-point segmentation} methods in a general nonparametric regression setting. The main contribution of this work is the adaptation property of these methods over a wide range of function classes, although they are designed for step functions. On the one hand, this includes the optimal convergence rates (up to log-factor) for step functions with bounded or even increasing to infinite number of jumps. On the other hand, for models beyond step functions, which are characterized by certain approximation spaces, we show the optimal rates (up to log-factor) as well. This includes bounded variation functions and (piecewise) H\"{o}lder functions of smoothness order $ 0 < \alpha \leq 1$. All results are formulated in terms of $L^p$-loss, $0 < p < \infty$, both almost surely and in expectation. In addition, we show that the convergence rates readily imply accuracy of feature detection, such as change-points, modes, troughs, etc. The practical performance is examined by various numerical simulations.

Dans le cadre du probleme de la regression parametrique de fonctions lisses, nous revisitons et analysons une approche variationnelle contrainte, que nous appelons l’estimateur de Nemirovski–Dantzig multi-echelle (MIND). Il peut etre vu comme une extension multi-echelle du selecteur de Dantzig (Ann. Statist. 35 (2009) 2313–2351), reposant sur des idees anterieures de Nemirovski (J. Comput. System Sci. 23 (1986) 1–11). MIND minimise une norme de Sobolev homogene sous la contrainte que la norme multi-resolution du terme residuel est bornee par un seuil universel. La contribution principale de cet article est d’obtenir des vitesses de convergence de MIND pour une norme $L^{q},1\leq q\leq \infty $, a la fois p.s. et en moyenne. A cette fin, nous introduisons la methode des conditions de sources approximees. Dans le cas d’un signal unidimensionnel, on peut exprimer ces dernieres en termes de proprietes d’approximation de B-splines. Une consequence remarquable est que MIND atteint presque des taux minimax optimaux simultanement pour un large ensemble d’espaces de Sobolev et Besov, ce qui montre une adaptabilite certaine. En plus de l’analyse asymptotique, nous etudions la performance de MIND pour un echantillon fini a l’aide de simulations numeriques. Un package MATLAB est disponible en ligne.

Recovering a function or high-dimensional parameter vector from indirect measurements is a central task in various scientific areas. Several methods for solving such inverse problems are well developed and well understood. Recently, novel algorithms using deep learning and neural networks for inverse problems appeared. While still in their infancy, these techniques show astonishing performance for applications like low-dose CT or various sparse data problems. However, there are few theoretical results for deep learning in inverse problems. In this paper, we establish a complete convergence analysis for the proposed NETT (network Tikhonov) approach to inverse problems. NETT considers nearly data-consistent solutions having small value of a regularizer defined by a trained neural network. We derive well-posedness results and quantitative error estimates, and propose a possible strategy for training the regularizer. Our theoretical results and framework are different from any previous work using neural networks for solving inverse problems. A possible data driven regularizer is proposed. Numerical results are presented for a tomographic sparse data problem, which demonstrate good performance of NETT even for unknowns of different type from the training data. To derive the convergence and convergence rates results we introduce a new framework based on the absolute Bregman distance generalizing the standard Bregman distance from the convex to the non-convex case.

Computing the inverse of pentadiagonal Toeplitz matrices has numerous applications. However, in some cases, we only need to obtain the diagonal of the inverse rather than the whole inverse. This paper presents a new approach to compute the diagonal of the pentadiagonal Toeplitz matrix inverse by Matrix Mobius Transformations. The proposed algorithm meticulously handles two cases of matrices. We have shown that the complexity of the proposed algorithm is O(n) for matrices of sizes n x n. We investigate the efficiency of the proposed algorithm via numerical simulations and make Matlab codes publicly available.

We introduce a new iterative regularization method for solving inverse problems that can be written as systems of linear or nonlinear equations in Hilbert spaces. The proposed averaged Kaczmarz (AVEK) method can be seen as a hybrid method between the Landweber and the Kaczmarz methods. As the Kaczmarz method, the proposed method only requires evaluation of one direct and one adjoint subproblem per iterative update. On the other hand, similar to the Landweber iteration, it uses an average over previous auxiliary iterates which increases stability. We present a convergence analysis of the AVEK iteration. Further, detailed numerical studies are presented for a tomographic image reconstruction problem, namely the limited data problem in photoacoustic tomography. Thereby, the AVEK is compared with other iterative regularization methods including standard Landweber and Kaczmarz iterations, as well as recently proposed accelerated versions based on error minimizing relaxation strategies.

Purpose In real-time MRI serial images are generally reconstructed from highly undersampled datasets as the iterative solutions of an inverse problem. While practical realizations based on regularized nonlinear inversion (NLINV) have hitherto been surprisingly successful, strong assumptions about the continuity of image features may affect the temporal fidelity of the estimated reconstructions. Theory and Methods The proposed method for real-time image reconstruction integrates the deformations between nearby frames into the data consistency term of the inverse problem. The aggregated motion estimation (AME) is not required to be affine or rigid and does not need additional measurements. Moreover, it handles multi-channel MRI data by simultaneously determining the image and its coil sensitivity profiles in a nonlinear formulation which also adapts to non-Cartesian (e.g., radial) sampling schemes. The new method was evaluated for real-time MRI studies using highly undersampled radial gradient-echo sequences. Results AME reconstructions for a motion phantom with controlled speed as well as for measurements of human heart and tongue movements demonstrate improved temporal fidelity and reduced residual undersampling artifacts when compared with NLINV reconstructions without motion estimation. Conclusion Nonlinear inverse reconstructions with aggregated motion estimation offer improved image quality and temporal acuity for visualizing rapid dynamic processes by real-time MRI. Magn Reson Med 72:1039–1048, 2014. © 2013 Wiley Periodicals, Inc.

The histogram is widely used as a simple, exploratory display of data, but it is usually not clear how to choose the number and size of bins. We construct a confidence set of distribution functions that optimally address the two main tasks of the histogram: estimating probabilities and detecting features such as increases and modes in the distribution. We define the essential histogram as the histogram in the confidence set with the fewest bins. Thus the essential histogram is the simplest visualization of the data that optimally achieves the main tasks of the histogram. The only assumption we make is that the data are independent and identically distributed. We provide a fast algorithm for the essential histogram, and illustrate our methodology with examples. An R-package is available on CRAN., Extension to discrete data is included. A R-package "essHist" is available from https://CRAN.R-project.org/package=essHist

Designing fast singular value decomposition (SVD) is significantly interesting in applications. The random direct SVD (RSVD) has provided a fast scheme to compute the well-approximate SVD by unilateral randomized sampling. In this paper, we present an efficient random algorithm in a bilateral sampling way. We also prove that the proposed algorithms can be bounded well and have less computational complexity compared to RSVD when the objective matrix is approximately square. Numerical experiments on graph Laplacian and Hilbert matrix demonstrate the efficiency and stability of the proposed methods.

Polynomials are widely used in scientific computing and engineering. In this paper, we present an accurate and fast compensated algorithm to evaluate bivariate polynomials with floating-point coefficients. This algorithm is applying error free transformations to the bivariate Horner scheme and sum the final decomposition accurately. We also prove the proposed algorithm's accuracy with forward error analysis that the accuracy of the computed result is similar to the result computed by the bivariate Horner scheme in twice the working precision. Numerical experiments illustrate the behavior and it has higher efficiency than the bivariate Horner scheme implemented in double-double library.

Fast multiple change-point segmentation methods, which additionally provide faithful statistical statements on the number, locations and sizes of the segments, have recently received great attention. In this paper, we propose a multiscale segmentation method, FDRSeg, which controls the false discovery rate (FDR) in the sense that the number of false jumps is bounded linearly by the number of true jumps. In this way, it adapts the detection power to the number of true jumps. We prove a non-asymptotic upper bound for its FDR in a Gaussian setting, which allows to calibrate the only parameter of FDRSeg properly. Change-point locations, as well as the signal, are shown to be estimated in a uniform sense at optimal minimax convergence rates up to a log-factor. The latter is w.r.t. $L^p$-risk, $p \ge 1$, over classes of step functions with bounded jump sizes and either bounded, or possibly increasing, number of change-points. FDRSeg can be efficiently computed by an accelerated dynamic program; its computational complexity is shown to be linear in the number of observations when there are many change-points. The performance of the proposed method is examined by comparisons with some state of the art methods on both simulated and real datasets. An R-package is available online., An R-package "FDRSeg" is available at http://www.stochastik.math.uni-goettingen.de/fdrs

This paper presents a compensated algorithm to accurately evaluate a polynomial expressed in Chebyshev basis of the first and second kind with floating-point coefficients. The principle is to apply error-free transformations to improve the traditional Clenshaw algorithm. The new algorithm is as accurate as the Clenshaw algorithm performed in twice the working precision. Forward error analysis and numerical experiments illustrate the accuracy and properties of the proposed algorithm.

Recent years have witnessed great success of manifold learning methods in understanding the structure of multidimensional patterns. However, most of these methods operate in a batch mode and cannot be effectively applied when data are collected sequentially. In this paper, we propose a general incremental learning framework, capable of dealing with one or more new samples each time, for the so-called spectral embedding methods. In the proposed framework, the incremental dimensionality reduction problem reduces to an incremental eigen-problem of matrices. Furthermore, we present, using this framework as a tool, an incremental version of Hessian eigenmaps, the IHLLE method. Finally, we show several experimental results on both synthetic and real world datasets, demonstrating the efficiency and accuracy of the proposed algorithm.

This paper presents a second-order and two-scale computation method of initial-boundary value problem for heat conduction equation with rapidly oscillating coefficients. This method can be applied to heat conduction analysis of composite materials. The main advantages of this method are that it can greatly save computer memory and CPU time, and it possesses good precision at the same time. Finally numerical results show that the presented method is effective and reliable.

Tensor-driven anisotropic diffusion and regularisation have been successfully applied to a wide range of image processing and computer vision tasks such as denoising, inpainting, and optical flow. Empirically it has been shown that anisotropic models with a diffusion tensor perform better than their isotropic counterparts with a scalar-valued diffusivity function. However, the reason for this superior performance is not well understood so far. Moreover, the specific modelling of the anisotropy has been carried out in a purely heuristic way. The goal of our paper is to address these problems. To this end, we use the statistics of natural images to derive a unifying framework for eight isotropic and anisotropic diffusion filters that have a corresponding variational formulation. In contrast to previous statistical models, we systematically investigate structure-adaptive statistics by analysing the eigenvalues of the structure tensor. With our findings, we justify existing successful models and assess the relationship between accurate statistical modelling and performance in the context of image denoising.