Scientists and clinicians interested in calibration-less parallel imaging techniques.Calibration-less parallel imaging methods [1-5] have been proposed to reduce the scan time needed to acquire the auto-calibration data in conventional calibration-based parallel imaging methods (e.g., GRAPPA [6] and SPIRiT [7]). Most of the existing calibration-less methods are based on low-rankness of matrices formed by the acquired k-space data. In this work, we propose a new formulation for calibration-free reconstruction using MultichAnneL Blind dEConvolution, named MALBEC. The method also formulates the reconstruction problem as a low rank matrix recovery problem, but the constraint is stronger than the existing SAKE method. Therefore, the proposed MALBEC is expected to provide better reconstructions. We formulate the reconstruction problem as a blind deconvolution problem in k-space, as illustrated in Figure 1. The acquired k-space data from K channels are denoted as $$${y_1}\left[{p,q}\right]$$$,$$${y_2}\left[{p,q}\right]$$$,...,$$${y_K}\left[{p,q}\right]$$$, which are the convolution between the k-space data $$$s[m,n]$$$ of the desired unknown image and the k-space of the coil sensitivities for each of the K different channels denoted as $$${h_1}\left[{m,n}\right]$$$,$$${h_2}\left[{m,n}\right]$$$,...,$$${h_k}\left[{m,n}\right]$$$. Mathematically, the convolution is given as $${y_k}\left[{p,q}\right]=\mathop\sum\limits_{m = 0}^{M-1}\mathop\sum\limits_{n=0}^{N-1}{h_k}\left[{m,n}\right]s\left[{p-m,q-n}\right].\quad(1)$$ In calibration-less parallel imaging, we wish to recover the full k-space data $$$s[m,n]$$$ from the undersampled $$${y_k}\left[{p,q}\right]$$$,$$$k=1,2,…,K$$$, without knowledge about $$$s[m,n]$$$ or $$${h_k}[m,n]$$$. Apparently, the problem is ill-posed with non-unique solutions. Since the coil sensitivities vary smoothly in image domain, to solve the ill-posed problem, we first assume their k-space values to have significant values only within a small window of size $$$M \times N$$$, which is much smaller than the size of $$$s$$$. As in Ref. [8], we then define a large rank-1 matrix $$${{\bf{X}}_0}={\bf{s}}[{\bf{h}}_1^{\rm{T}}{\bf{h}}_2^{\rm{T}}\cdots{\bf{h}}_K^{\rm{T}}]$$$ where $$$s$$$ is the vectorized $$$s[-m,-n]$$$ and $$${{\bf{h}}_k}$$$ vectorized. It is seen that $$${{\bf{X}}_0}$$$ contains all pairs of variables appearing in the sum in Eq. (1). As illustrated in Figure 2 (using a 1D convolution example), each acquired k-space data point $$${y_k}\left[{p,q}\right]$$$ is now a linear combination of the entries in $$${{\bf{X}}_0}$$$ along the skew direction. We define a linear operation $$${\cal A}(\cdot)$$$ which sums over skew diagonals of the submatrices $$${\bf{sh}}_k^{\rm{T}}$$$ in $$${{\bf{X}}_0}$$$ to generate $$${y_k}\left[ {p,q} \right]$$$, and define as the undersampling operation. Concatenating the undersampled observations from all channels into a single vector $$${{\bf{y}}_{\rm{\Omega }}}$$$, we have a linear system of equations: $${{\bf{y}}_{\rm{\Omega }}} = {\rm{\Omega }}{\cal A}({{\bf{X}}_0}). \quad(2)$$ The problem of reconstructing $$$s$$$ now becomes a matrix recovery problem of $$${{\bf{X}}_0}$$$. Because we know that the unknown matrix $$${{\bf{X}}_0}$$$ has rank one, to recover $$${{\bf{X}}_0}$$$ from $$${{\bf{y}}_{\rm{\Omega }}}$$$, we want to find a matrix that satisfies Eq. (2) such that $$$rank({{\bf{X}}_0}) = 1$$$. However, the problem is difficult to solve directly. We therefore solve $$${\min_{\bf{X}}}{\left\| {\bf{X}} \right\|_{\rm{*}}}{\rm{}}$$$ subject to $$${{\bf{y}}_{\rm{\Omega }}} = {\rm{\Omega }}{\cal A}({\bf{X}}) \quad(3)$$$ instead, where $$${\left\| {\bf{X}}\right\|_*}$$$ is the nuclear norm (sum of the singular values) of $$${\bf{X}}$$$. It has been proved that under certain conditions (including randomness in $$${\rm{\Omega }}$$$, Eq. (3) can (with high probability) recover exactly the rank-1 matrix. Conjugate gradient method was used to solve $$$\left( {{{\cal A}^{ - 1}}{{\rm{\Omega }}^{ - 1}}{\rm{\Omega }}{\cal A}+{{\rm{\lambda}}^{-1}}{\rm{{\rm I}}}} \right){\bf{X}}={{\rm{\lambda }}^{-1}}{\bf{G}}+{{\cal A}^{-1}}{{\rm{\Omega}}^{-1}}{{\bf{y}}_{\rm{\Omega}}}$$$, where $$${\bf{G}} = {{\bf{u}}_1}{{\bf{v}}_1}$$$, $$${{\bf{u}}_1}$$$ and $$${{\bf{v}}_1}$$$ are the first left and right singular vectors of matrix $$${\bf{X}}$$$. Once we have $$${\bf{X}}$$$, we can find the full k-space data at all channels using $$${\bf{y}} = {\cal A}({\bf{X}})$$$. The desired image can be obtained as the sum-of squares of the images from all channels. To evaluate the performance of the proposed method, we use the human brain dataset in Ref. [1] (3D SPGR sequence on a 1.5T MRI scanner (GE, Waukesha, WI) with an 8-channel head coil, TE=8 ms, TR=17.6 ms, and flip angle=20°. FOV=20cm×20cm×20cm, matrix size of 200×200×200). A single slice was selected from this data set and was used throughout the experiments. The data were fully acquired and then retrospectively and randomly under-sampled (with sampling pattern in Ref. [9]) to simulate the accelerated acquisition with a reduction factor of 3. The proposed MALBEC is compared against another rank-constrained method, SAKE [1]. Figure 3 shows that the MALBEC reconstruction quality is better than that of SAKE, while the computational time is about the same. SAKE can be viewed as a relaxed version of the proposed method. Specifically, using the same notation, SAKE can be formulated as $$${\rm{mi}}{{\rm{n}}_{\bf{X}}}{\rm{rank}}((1+{\bf{P}}+\cdots+{{\bf{P}}^{MN}}){\bf{X}})$$$ subject to $$${{\bf{y}}_{\rm{\Omega}}} = {\rm{\Omega }}{\cal A}\left({\bf{X}}\right)$$$, where $$${\bf{P}}$$$ is a circulant shift matrix. It is seen that this constraint is weaker than that in Eq. (3).

A novel calibration-less parallel imaging technique is proposed. It is shown to be superior to SAKE in both theory and experiments.