Golub-Kahan Bidiagonalization Research Articles

Robust residual-guided iterative reconstruction for sparse-view CT in small animal imaging

Objective. We introduce a robust image reconstruction algorithm named residual-guided Golub–Kahan iterative reconstruction technique (RGIRT) designed for sparse-view computed tomography (CT), which aims at high-fidelity image reconstruction from a limited number of projection views. Approach. RGIRT utilizes an inner-outer dual iteration framework, with a flexible least square QR (FLSQR) algorithm implemented in the inner iteration and a restarted iterative scheme applied in the outer iteration. The inner FLSQR employs a flexible Golub–Kahan bidiagonalization method to reduce the size of the inverse problem, and a weighted generalized cross-validation method to adaptively estimate the regularization hyper-parameter. The inner iteration efficiently yields the intermediate reconstruction result, while the outer iteration minimizes the residual and refines the solution by using the result obtained from the inner iteration. Main results. The reconstruction performance of RGIRT is evaluated and compared to other reference methods (FBPConvNet, SART-TV, and FLSQR) using projection data from both numerical phantoms and real experimental Micro-CT data. The experimental findings, from testing various numbers of projection views and different noise levels, underscore the robustness of RGIRT. Meanwhile, theoretical analysis confirms the convergence of residual for our approach. Significance. We propose a robust iterative reconstruction algorithm for x-ray CT scans with sparse views, thereby shortening scanning time and mitigating excessive ionizing radiation exposure to small animals.

A CCBM-based generalized GKB iterative regularization algorithm for inverse Cauchy problems

This paper examines inverse Cauchy problems that are governed by a kind of elliptic partial differential equation. The inverse problems involve recovering the missing data on an inaccessible boundary from the measured data on an accessible boundary, which is severely ill-posed. By using the coupled complex boundary method (CCBM), which integrates both Dirichlet and Neumann data into a single Robin boundary condition, we reformulate the underlying problem into an operator equation. Based on this new formulation, we study the solution existence issue of the reduced problem with noisy data. A Golub–Kahan bidiagonalization (GKB) process together with Givens rotation is employed for iteratively solving the proposed operator equation. The regularizing property of the developed method, called CCBM-GKB, and its convergence rate results are proved under a posteriori stopping rule. Finally, a linear finite element method is used for the numerical realization of CCBM-GKB. Various numerical experiments demonstrate that CCBM-GKB is a kind of accelerated iterative regularization method, as it is much faster than the classic Landweber method.

Krylov Subspace Approach to Core Problems within Multilinear Approximation Problems: A Unifying Framework

Error contaminated linear approximation problems appear in a large variety of applications. The presence of redundant or irrelevant data complicates their solution. It was shown that such data can be removed by the core reduction yielding a minimally dimensioned subproblem called the core problem. Direct (SVD or Tucker decomposion-based) reduction has been introduced previously for problems with matrix models and vector, or matrix, or tensor observations; and also for problems with bilinear models. For the cases of vector and matrix observations a Krylov subspace method, the generalized Golub–Kahan bidiagonalization, can be used to extract the core problem. In this paper, we first unify previously studied variants of linear approximation problems under the general framework of a multilinear approximation problem. We show how the direct core reduction can be extended to it. Then we show that the generalized Golub–Kahan bidiagonalization yields the core problem for any multilinear approximation problem. This further allows one to prove various properties of core problems, in particular, we give upper bounds on the multiplicity of singular values of reduced matrices.

A Majorization-Minimization Golub-Kahan bidiagonalization method for <inline-formula><tex-math id="M1">$ \ell_{2}-\ell_{q} $</tex-math></inline-formula> mimimization with applications in image restorization

An image restorization problem is often modelled as a discrete ill-posed problem. In general, its solution, even if it exists, is very sensitive to the perturbation in the data. Regularization methods reduces the sensitivity by replacing this problem with a minimization problem with a fidelity term and $ \ell_{q} $ regularization term. In order to improve the sparsity of the solution, we only consider the case of $ 0<q\leq1 $ in this paper. This paper presents a majorization-minimization Golub-Kahan bidiagonalization algorithm to solve this kind of minimization problems. The solution subspace is extended by the Golub-Kahan bidiagonalization process. The restarted case is also considered. The regularization parameter is determined by using the discrepancy principle. Several examples in image restorization are shown for the proposed methods.

Error estimates for Golub–Kahan bidiagonalization with Tikhonov regularization for ill–posed operator equations

Linear ill-posed operator equations arise in various areas of science and engineering. The presence of errors in the operator and the data often makes the computation of an accurate approximate solution difficult. In this paper, we compute an approximate solution of an ill-posed operator equation by first determining an approximation of the operators of generally fairly small dimension by carrying out a few steps of a continuous version of the Golub–Kahan bidiagonalization process to the noisy operator. Then Tikhonov regularization is applied to the low-dimensional problem so obtained and the regularization parameter is determined by solving a low-dimensional nonlinear equation. The effect of the errors incurred in each step of the solution process is analyzed. Computed examples illustrate the theory presented.

Inexact inner–outer Golub–Kahan bidiagonalization method: A relaxation strategy

AbstractWe study an inexact inner–outer generalized Golub–Kahan algorithm for the solution of saddle‐point problems with a two‐times‐two block structure. In each outer iteration, an inner system has to be solved which in theory has to be done exactly. Whenever the system is getting large, an inner exact solver is, however, no longer efficient or even feasible and iterative methods must be used. We focus this article on a numerical study showing the influence of the accuracy of an inner iterative solution on the accuracy of the solution of the block system. Emphasis is further given on reducing the computational cost, which is defined as the total number of inner iterations. We develop relaxation techniques intended to dynamically change the inner tolerance for each outer iteration to further minimize the total number of inner iterations. We illustrate our findings on a Stokes problem and validate them on a mixed formulation of the Poisson problem.

Computation of generalized matrix functions with rational Krylov methods

We present a class of algorithms based on rational Krylov methods to compute the action of a generalized matrix function on a vector. These algorithms incorporate existing methods based on the Golub-Kahan bidiagonalization as a special case. By exploiting the quasiseparable structure of the projected matrices, we show that the basis vectors can be updated using a short recurrence, which can be seen as a generalization to the rational case of the Golub-Kahan bidiagonalization. We also prove error bounds that relate the error of these methods to uniform rational approximation. The effectiveness of the algorithms and the accuracy of the bounds is illustrated with numerical experiments.

Tensor Krylov subspace methods via the Einstein product with applications to image and video processing

In the present paper, we are interested in developing iterative Krylov subspace methods in tensor structure to solve a class of multilinear systems via the Einstein product. In particular, we develop tensor variants of the GMRES and Golub–Kahan bidiagonalization processes in tensor framework. We further consider the case that mentioned equation may be possibly corresponds to a discrete ill-posed multidimensional problem. Applications arising from color image and video restoration are included.

When does the Lanczos algorithm compute exactly?

In theory, the Lanczos algorithm generates an orthogonal basis of the corresponding Krylov subspace. However, in finite precision arithmetic the orthogonality and linear independence of the computed Lanczos vectors is usually lost quickly. In this paper we study a class of matrices and starting vectors having a special nonzero structure that guarantees exact computations of the Lanczos algorithm whenever floating point arithmetic satisfying the IEEE 754 standard is used. Analogous results are formulated also for an implementation of the conjugate gradient method called cgLanczos. This implementation then computes approximations that agree with their exact counterparts to a relative accuracy given by the machine precision and the condition number of the system matrix. The results are extended to the Arnoldi algorithm, the nonsymmetric Lanczos algorithm, the Golub-Kahan bidiagonalization, the block-Lanczos algorithm, and their counterparts for solving linear systems.

The tensor Golub–Kahan–Tikhonov method applied to the solution of ill‐posed problems with a t‐product structure

AbstractThis paper discusses an application of partial tensor Golub–Kahan bidiagonalization to the solution of large‐scale linear discrete ill‐posed problems based on the t‐product formalism for third‐order tensors proposed by Kilmer and Martin (M. E. Kilmer and C. D. Martin, Factorization strategies for third order tensors, Linear Algebra Appl., 435 (2011), pp. 641‐658). The solution methods presented first reduce a given (large‐scale) problem to a problem of small size by application of a few steps of tensor Golub–Kahan bidiagonalization and then regularize the reduced problem by Tikhonov's method. The regularization operator is a third‐order tensor, and the data may be represented by a matrix, that is, a tensor slice, or by a general third‐order tensor. A regularization parameter is determined by the discrepancy principle. This results in fully automatic solution methods that neither require a user to choose the number of bidiagonalization steps nor the regularization parameter. The methods presented extend available methods for the solution for linear discrete ill‐posed problems defined by a matrix operator to linear discrete ill‐posed problems defined by a third‐order tensor operator. An interlacing property of singular tubes for third‐order tensors is shown and applied. Several algorithms are presented. Computed examples illustrate the advantage of the tensor t‐product approach, in comparison with solution methods that are based on matricization of the tensor equation.

Accurate and fast matrix factorization for low-rank learning.

In this paper, we tackle two important problems in low-rank learning, which are partial singular value decomposition and numerical rank estimation of huge matrices. By using the concepts of Krylov subspaces such as Golub-Kahan bidiagonalization (GK-bidiagonalization) as well as Ritz vectors, we propose two methods for solving these problems in a fast and accurate way. Our experiments show the advantages of the proposed methods compared to the traditional and randomized singular value decomposition methods. The proposed methods are appropriate for applications involving huge matrices where the accuracy of the desired singular values and also all of their corresponding singular vectors are essential. As a real application, we evaluate the performance of our methods on the problem of Riemannian similarity learning between two different image datasets of MNIST and USPS.

On tensor GMRES and Golub-Kahan methods via the T-product for color image processing

The present paper is concerned with developing tensor iterative Krylov subspace methods to solve large multi-linear tensor equations. We use the T-product for two tensors to define tensor tubal global Arnoldi and tensor tubal global Golub-Kahan bidiagonalization algorithms. Furthermore, we illustrate how tensor-based global approaches can be exploited to solve ill-posed problems arising from recovering blurry multichannel (color) images and videos, using the so-called Tikhonov regularization technique, to provide computable approximate regularized solutions. We also review a generalized cross-validation and discrepancy principle type of criterion for the selection of the regularization parameter in the Tikhonov regularization. Applications to image sequence processing are given to demonstrate the efficiency of the algorithms.

Golub\u2013Kahan vs. Monte Carlo: a comparison of bidiagonlization and a randomized SVD method for the solution of linear discrete ill-posed problems

Randomized methods can be competitive for the solution of problems with a large matrix of low rank. They also have been applied successfully to the solution of large-scale linear discrete ill-posed problems by Tikhonov regularization (Xiang and Zou in Inverse Probl 29:085008, 2013). This entails the computation of an approximation of a partial singular value decomposition of a large matrix A that is of numerical low rank. The present paper compares a randomized method to a Krylov subspace method based on Golub–Kahan bidiagonalization with respect to accuracy and computing time and discusses characteristics of linear discrete ill-posed problems that make them well suited for solution by a randomized method.

On the block Lanczos and block Golub–Kahan reduction methods applied to discrete ill‐posed problems

AbstractThe reduction of a large‐scale symmetric linear discrete ill‐posed problem with multiple right‐hand sides to a smaller problem with a symmetric block tridiagonal matrix can easily be carried out by the application of a small number of steps of the symmetric block Lanczos method. We show that the subdiagonal blocks of the reduced problem converge to zero fairly rapidly with increasing block number. This quick convergence indicates that there is little advantage in expressing the solutions of discrete ill‐posed problems in terms of eigenvectors of the coefficient matrix when compared with using a basis of block Lanczos vectors, which are simpler and cheaper to compute. Similarly, for nonsymmetric linear discrete ill‐posed problems with multiple right‐hand sides, we show that the solution subspace defined by a few steps of the block Golub–Kahan bidiagonalization method usually can be applied instead of the solution subspace determined by the singular value decomposition of the coefficient matrix without significant, if any, reduction of the quality of the computed solution.

Hybrid projection methods for large-scale inverse problems with mixed Gaussian priors

When solving ill-posed inverse problems, a good choice of the prior is critical for the computation of a reasonable solution. A common approach is to include a Gaussian prior, which is defined by a mean vector and a symmetric and positive definite covariance matrix, and to use iterative projection methods to solve the corresponding regularized problem. However, a main challenge for many of these iterative methods is that the prior covariance matrix must be known and fixed (up to a constant) before starting the solution process. In this paper, we develop hybrid projection methods for inverse problems with mixed Gaussian priors where the prior covariance matrix is a convex combination of matrices and the mixing parameter and the regularization parameter do not need to be known in advance. Such scenarios may arise when data is used to generate a sample prior covariance matrix (e.g., in data assimilation) or when different priors are needed to capture different qualities of the solution. The proposed hybrid methods are based on a mixed Golub–Kahan process, which is an extension of the generalized Golub–Kahan bidiagonalization, and a distinctive feature of the proposed approach is that both the regularization parameter and the weighting parameter for the covariance matrix can be estimated automatically during the iterative process. Furthermore, for problems where training data are available, various data-driven covariance matrices (including those based on learned covariance kernels) can be easily incorporated. Numerical examples from tomographic reconstruction demonstrate the potential for these methods.

Two-Level Iterative Methods for Solving the Saddle Point Problems

Iterative processes in the Krylov subspaces for solving large ill conditioned saddle-type SLAEs with sparse matrices arising in finite difference, finite volume, and finite element approximations of multidimensional boundary value problems with complex geometric and functional properties of the initial data, characteristic of many relevant applications are studied. Combined two-level iterative algorithms using efficient Chebyshev acceleration and variational the conjugate directions methods, as well as the Golub-Kahan bi-diagonalization algorithms in the Krylov subspaces are considered. Examples of two-dimensional and three-dimensional filtration problems are used to study the resource consumption and computational performance of the proposed algorithms, as well as their scalable parallization on the multiprocessor systems with distributed and hierarchical shared memory.

Application of an iterative Golub-Kahan algorithm to structural mechanics problems with multi-point constraints

Kinematic relationships between degrees of freedom, also named multi-point constraints, are frequently used in structural mechanics. In this paper, the Craig variant of the Golub-Kahan bidiagonalization algorithm is used as an iterative method to solve the arising linear system with a saddle point structure. The condition number of the preconditioned operator is shown to be close to unity and independent of the mesh size. This property is proved theoretically and illustrated on a sequence of test problems of increasing complexity, including concrete structures enforced with pretension cables and the coupled finite element model of a reactor containment building. The Golub-Kahan algorithm converges in only a small number of steps for all considered test problems and discretization sizes. Furthermore, it is robust in practical cases that are otherwise considered to be difficult for iterative solvers.

A fast methodology for large-scale focusing inversion of gravity and magnetic data using the structured model matrix and the 2-D fast Fourier transform

SUMMARYWe discuss the focusing inversion of potential field data for the recovery of sparse subsurface structures from surface measurement data on a uniform grid. For the uniform grid, the model sensitivity matrices have a block Toeplitz Toeplitz block structure for each block of columns related to a fixed depth layer of the subsurface. Then, all forward operations with the sensitivity matrix, or its transpose, are performed using the 2-D fast Fourier transform. Simulations are provided to show that the implementation of the focusing inversion algorithm using the fast Fourier transform is efficient, and that the algorithm can be realized on standard desktop computers with sufficient memory for storage of volumes up to size n ≈ 106. The linear systems of equations arising in the focusing inversion algorithm are solved using either Golub–Kahan bidiagonalization or randomized singular value decomposition algorithms. These two algorithms are contrasted for their efficiency when used to solve large-scale problems with respect to the sizes of the projected subspaces adopted for the solutions of the linear systems. The results confirm earlier studies that the randomized algorithms are to be preferred for the inversion of gravity data, and for data sets of size m it is sufficient to use projected spaces of size approximately m/8. For the inversion of magnetic data sets, we show that it is more efficient to use the Golub–Kahan bidiagonalization, and that it is again sufficient to use projected spaces of size approximately m/8. Simulations support the presented conclusions and are verified for the inversion of a magnetic data set obtained over the Wuskwatim Lake region in Manitoba, Canada.

Regularization properties of LSQR for linear discrete ill-posed problems in the multiple singular value case and best, near best and general low rank approximations**The work was supported by the National Science Foundation of China (No. 11771249).

For the large-scale linear discrete ill-posed problem min‖Ax − b‖ or Ax = b with b contaminated by white noise, the Golub–Kahan bidiagonalization based LSQR method and its mathematically equivalent CGLS, the conjugate gradient (CG) method applied to ATAx = ATb, are most commonly used. They have intrinsic regularizing effects, where the iteration number k plays the role of regularization parameter. The long-standing fundamental question is: Can LSQR and CGLS find two-norm filtering best possible regularized solutions? The author has given definitive answers to this question for severely and moderately ill-posed problems when the singular values of A are simple. This paper extends the results to the multiple singular value case, and studies the approximation accuracy of Krylov subspaces, the quality of low rank approximations generated by Golub–Kahan bidiagonalization and the convergence properties of Ritz values. For the two kinds of problems, we prove that LSQR finds two-norm filtering best possible regularized solutions at semi-convergence. Particularly, we consider some important and untouched issues on best, near best and general rank k approximations to A for the ill-posed problems with the singular values with α > 0, and the relationships between them and their nonzero singular values. Numerical experiments confirm our theory. The results on general rank k approximations and the properties of their nonzero singular values apply to several Krylov solvers, including LSQR, CGME, MINRES, MR-II, GMRES and RRGMRES.

A right preconditioner for the LSMR method

The LSMR (Least Squares Minimal Residual) method is an absorbing solver that can solve linear system Ax = b and least squares problem min ||Ax = b|| where A is a sparse and large matrix. This method is based on the Golub-Kahan bidiagonalization process and sometimes it may converge slowly like other methods. I n order to prevent this event, a right preconditioner for LSMR method is presented to solve large and sparse linear system which used for LSQR (Least Squares with QR factorization) method before. Numerical examples and comparing the preconditioned LSMR method to unpreconditioned LSMR method would show the effectiveness of the preconditioner. I t is obtained from this paper that PLSMR (Preconditioned LSMR) method has a better performance in reducing the number of iterations and relative residual norm in comparing with the original LSMR method.