1

Convolutional Recurrent Neural Networks for
Dynamic MR Image Reconstruction

arXiv:1712.01751v2 [] 26 Jan 2018

Chen Qin*† , Jo Schlemper*, Jose Caballero, Anthony N. Price, Joseph V. Hajnal and Daniel
Rueckert Fellow, IEEE

In order to accelerate MRI acquisition, most approaches
consider undersampling the data in k-space (frequency domain). Due to the violation of the Nyquist sampling theorem,
undersampling introduces aliasing artefacts in the image
domain. Images can be subsequently reconstructed by solving
an optimisation problem that regularises the solution with
assumptions on the underlying data, such as smoothness,
sparsity or, for the case of dynamic imaging, spatio-temporal
redundancy. Past literature has shown that exploiting spatiotemporal redundancy can greatly improve image reconstruction
quality compared to compressed sensing (CS) based single
frame reconstruction methods [1], [2]. However, the limitations
of these optimisation based approaches are the following: firstly,
it is non-trivial to propose an optimisation function without
introducing significant bias on the considered data. In addition,
the manual adjustments of hyperparameters are nontrivial.
Secondly, these optimisation problems often involve highly
nonlinear, non-convex terms. As a consequence, the majority
of approaches resort on iterative algorithms to reconstruct the
images, in which often neither attaining the global minimum
nor convergence to a solution is guaranteed. Furthermore,
the reconstruction speeds of these methods are often slow.
Proposing a robust iterative algorithm is still an active area of
Index Terms—Recurrent neural network, convolutional neural research.
network, dynamic magnetic resonance imaging, cardiac image
In comparison, deep learning methods are gaining popularity
reconstruction
for their accuracy and efficiency. Unlike traditional approaches,
the prior information and regularisation are learnt implicitly
from data, allowing the reconstruction to look more natural.
I. I NTRODUCTION
However, the limitation is that so far only a handful of
AGNETIC Resonance Imaging (MRI) is a non-invasive approaches exist [3], [4] for dynamic reconstruction. Hence, the
imaging technique which offers excellent spatial reso- applicability of deep learning models to this problem is yet to be
lution and soft tissue contrast and is widely used for clinical fully explored. In particular, a core question is how to optimally
diagnosis and research. Dynamic MRI attempts to reveal both exploit spatio-temporal redundancy. In addition, most deep
spatial and temporal profiles of the underlying anatomy, which learning methods do not exploit domain-specific knowledge,
has a variety of applications such as cardiovascular imaging which potentially enables the networks to efficiently learn data
and perfusion imaging. However, the acquisition speed is representation by regulating the mechanics of network layers,
fundamentally limited due to both hardware and physiological hence boosting their performances.
constraints as well as the requirement to satisfy the Nyquist
In this work, we propose a novel convolutional recurrent
sampling rate. Long acquisition times are not only a burden neural network (CRNN) method to reconstruct high quality
for patients but also make MRI susceptible to motion artefacts. dynamic MR image sequences from undersampled data, termed
CRNN-MRI, which aims to tackle the aforementioned prob† Corresponding author: Chen Qin. (Email address: c.qin15@imperial.ac.uk)
lems of both traditional and deep learning methods. Firstly,
*These authors contributed equally to this work.
we
formulate a general optimisation problem for solving
C. Qin, J. Schlemper, J. Caballero and D. Rueckert are with the Biomedical
Image Analysis Group, Department of Computing, Imperial College London, accelerated dynamic MRI based on variable splitting and
SW7 2AZ London, UK.
alternate minimisation. We then show how this algorithm
J. V. Hajnal, and A. N. Price are with with the Division of Imaging Sciences
and Biomedical Engineering Department, King’s College London, St. Thomas’ can be seen as a network architecture. In particular, the
Hospital, SE1 7EH London, U.K.
proposed method consists of a CRNN block which acts as the

Abstract—Accelerating the data acquisition of dynamic magnetic resonance imaging (MRI) leads to a challenging ill-posed
inverse problem, which has received great interest from both
the signal processing and machine learning communities over
the last decades. The key ingredient to the problem is how
to exploit the temporal correlations of the MR sequence to
resolve aliasing artefacts. Traditionally, such observation led
to a formulation of a non-convex optimisation problem, which
was solved using iterative algorithms. Recently, however, deep
learning based-approaches have gained significant popularity
due to their ability to solve general inverse problems. In this
work, we propose a unique, novel convolutional recurrent neural
network (CRNN) architecture which reconstructs high quality
cardiac MR images from highly undersampled k-space data by
jointly exploiting the dependencies of the temporal sequences
as well as the iterative nature of the traditional optimisation
algorithms. In particular, the proposed architecture embeds
the structure of the traditional iterative algorithms, efficiently
modelling the recurrence of the iterative reconstruction stages
by using recurrent hidden connections over such iterations. In
addition, spatio-temporal dependencies are simultaneously learnt
by exploiting bidirectional recurrent hidden connections across
time sequences. The proposed method is able to learn both the
temporal dependency and the iterative reconstruction process
effectively with only a very small number of parameters, while
outperforming current MR reconstruction methods in terms of
computational complexity, reconstruction accuracy and speed.

M

2

proximal operator and a data consistency layer corresponding Furthermore, these methods are not able to exploit the prior
to the classical data fidelity term. In addition, the CRNN knowledge that can be learnt from the vast number of MRI
block employs recurrent connections across each iteration step, exams routinely performed, which should be helpful to further
allowing reconstruction information to be shared across the guide the reconstruction process.
multiple iterations of the process. Secondly, we incorporate
Recently, deep learning-based MR reconstruction has gained
bidirectional convolutional recurrent units evolving over time popularity due to its promising results for solving inverse and
to exploit the temporal dependency of the dynamic sequences compressed sensing problems. In particular, two paradigms
and effectively propagate the contextual information across have emerged: the first class of approaches proposes to use
time frames of the input. As a consequence, the unique CRNN convolutional neural networks (CNNs) to learn an end-to-end
architecture jointly learns representations in a recurrent fashion mapping, where architectures such as SRCNN [10] or U-net
evolving over both time sequences as well as iterations of the [11] are often chosen for MR image reconstruction [12], [13],
reconstruction process, effectively combining the benefits of [14], [15]. The second class of approaches attempts to make
traditional iterative methods and deep learning.
each stage of iterative optimisation learnable by unrolling the
To the best of our knowledge, this is the first work applying end-to-end pipeline into a deep network [9], [16], [17], [18],
RNNs for dynamic MRI reconstruction. The contributions of [19]. For instance, Hammernik et al. [9] introduced a trainable
this work are the following: Firstly, we view the optimisation formulation for accelerated parallel imaging (PI) based MRI
problem of dynamic data as a recurrent network and describe reconstruction termed variational network, which embedded a
a novel CRNN architecture which simultaneously incorporates CS concept within a deep learning approach. ADMM-Net [17]
the recurrent relationship of data over time and iterations. was proposed by reformulating an alternating direction method
Secondly, we demonstrate that the proposed method shows of multipliers (ADMM) algorithm to a deep network, where
promising results and improves upon the current state-of-the-art each stage of the architecture corresponds to an iteration in
dynamic MR reconstruction methods both in reconstruction the ADMM algorithm. More recently, Schlemper et al. [18]
accuracy and speed. Finally, we compare our architecture to 3D proposed a cascade network which simulated the iterative
CNN which does not impose the recurrent structure. We show reconstruction of dictionary learning-based methods and were
that the proposed method outperforms the CNN at different later extended for dynamic MR reconstructions [3]. Most
undersampling rates and speed, while requiring significantly approaches so far have focused on 2D images, whereas only a
few approaches exist for dynamic MR reconstruction [3], [4].
fewer parameters.
While they show promising results, the optimal architecture,
training scheme and configuration spaces are yet to be fully
II. R ELATED W ORK
explored.
One of the main challenges associated with recovering an
More recently, several ideas were proposed, which integrate
uncorrupted image is that both the undersampling strategy deep neural networks with model-based optimisation methods
and a-priori knowledge of appropriate properties of the image to solve inverse problems [20], [21]. In contrast to these papers,
need to be taken into account. Methods like k-t BLAST and which proposed to deal with a fidelity term and a regularisation
k-t SENSE [5] take advantage of a-priori information about term separately, we integrate the two terms in a single deep
the x-f support obtained from the training data set in order network, so that the network can be trained end-to-end. In
to prune a reconstruction to optimally reduce aliasing. An addition, as we will show, we integrate a hidden connection
alternative popular approach is to exploit temporal redundancy over the optimisation iteration to enable the information used
to unravel from the aliasing by using CS approaches [1], [6] or for the reconstruction at each iteration to be shared across all
CS combined with low-rank approaches [2], [7]. The class of stages of the reconstruction process, aiming for an iterative
methods which employ CS to the MRI reconstruction is termed algorithm that can fully benefit from information extracted at
as CS-MRI [8]. They assume that the image to be reconstructed all processing stages. As to the nature of the proposed RNN
has a sparse representation in a certain transform domain, and units, previous work involving RNNs only updated the hidden
they need to balance sparsity in the transform domain against state of the recurrent connection with a fixed input [22], [23],
consistency with the acquired undersampled k-space data. For [24], while the proposed architecture progressively updates
instance, an example of successful methods enforcing sparsity the input as the optimisation iteration increases. In addition,
in x-f domain is k-t FOCUSS [1]. A low rank and sparse previous work only modelled the recurrence of iteration or time
reconstruction scheme (k-t SLR) [2] introduces non-convex [25] exclusively, whereas the proposed method jointly exploits
spectral norms and uses a spatio-temporal total variation norm both dimensions, yielding a unique architecture suitable for
in recovering the dynamic signal matrix. Dictionary learning the dynamic reconstruction problem.
approaches were also proposed to train an over-complete basis
of atoms to optimally sparsify spatio-temporal data [6]. These
III. C ONVOLUTIONAL R ECURRENT N EURAL N ETWORK
methods offer great potential for accelerated imaging, however,
FOR MRI RECONSTRUCTION
they often impose strong assumptions on the underlying data,
requiring nontrivial manual adjustments of hyperparameters A. Problem Formulation
depending on the application. In addition, it has been observed
Let x ∈ CD denote a sequence of complex-valued MR
that these methods tend to result in blocky [9] and unnatural images to be reconstructed, represented as a vector with
reconstructions, and their reconstruction speed is often slow. D = Dx Dy T , and let y ∈ CM (M << D) represent the

3

undersampled k-space measurements, where Dx and Dy are
width and height of the frame respectively and T stands
for the number of frames. Our problem is to reconstruct x
from y, which is commonly formulated as an unconstrained
optimisation problem of the form:
argmin R(x) + λky − Fu xk22
x

(1)

Here Fu is an undersampling Fourier encoding matrix, R
expresses regularisation terms on x and λ allows the adjustment
of data fidelity based on the noise level of the acquired
measurements y. For CS and low-rank based approaches, the
regularisation terms R often employed are `0 or `1 norms in
the sparsifying domain of x as well as the rank or nuclear
norm of x respectively. In general, Eq. 1 is a non-convex
function and hence, the variable splitting technique is usually
adopted to decouple the fidelity term and the regularisation
term. By introducing an auxiliary variable z that is constrained
to be equal to x, Eq. 1 can be reformulated to minimize the
following cost function via the penalty method:
argmin R(z) + λky − Fu xk22 + µkx − zk22

(2)

individually parameterise each step. In our formulation, we
model each optimisation stage (i) as a learnt, recurrent, forward
encoding step fi (x(i−1) , z(i−1) ; θ, y, λ, Ω). The difference
is that now we use one model which performs proximal
operator, however, it also allows itself to propagate information
across iteration, making it adaptable for the changes across the
optimisation steps. The detail will be discussed in the following
section. The different strategies are illustrated in Fig 1.
B. CRNN for MRI reconstruction
RNN is a class of neural networks that makes use of
sequential information to process sequences of inputs. They
maintain an internal state of the network acting as a "memory",
which allows RNNs to naturally lend themselves to the
processing of sequential data. Inspired by iterative optimisation
schemes of Eq. 3, we propose a novel convolutional RNN
(CRNN) network. In the most general scope, our neural
encoding model is defined as follows,
xrec = fN (fN −1 (· · · (f1 (xu )))),

(5)

x,z

where µ is a penalty parameter. By applying alternate minimisation over x and z, Eq. 2 can be solved via the following
iterative procedures:
z(i+1) = argmin R(z) + µkx(i) − zk22

(3a)

x(i+1) = argmin λky − Fu xk22 + µkx − z(i+1) k22

(3b)

z

in which xrec denotes the prediction of the network, xu is the
sequence of undersampled images with length T and also the
input of the network, fi (xu ; θ, λ, Ω) is the network function
for each iteration of optimisation step, and N is the number
of iterations. We can compactly represent a single iteration fi
of our network as follows:

x

where x(0) = xu = FH
u y is the zero-filled reconstruction taken
as an initialisation and z can be seen as an intermediate state
of the optimisation process. For MRI reconstruction, Eq. 3b is
often regarded as a data consistency (DC) step where we can
obtain the following closed-form solution [18]:

(i−1)
x(i)
+ CRNN(x(i−1)
rnn = xrec
rec ),

x(i)
rec

=

DC(x(i)
rnn ; y, λ0 , Ω),

(6a)
(6b)

where CRNN is a learnable block explained hereafter, DC
(i)
is the data consistency step treated as a network layer, xrec
is the progressive reconstruction of the undersampled image
λ0
H
(i)
H
(i)
(0)
(i)
x(i+1) =
xu at iteration i with xrec = xu , xrnn is the intermediate
( DC(z ; y, λ0 , Ω) = F ΛFz + 1+λ0 Fu y,
reconstruction image before the DC layer, and y is the
1
if k 6∈ Ω
Λkk =
1
acquired k-space samples. Note that the variables xrec , xrnn are
if
k
∈
Ω
1+λ0
(4) analogous to x, z in Eq. 3 respectively. Here, we use CRNN
in which F is the full Fourier encoding matrix (a discrete to encode the update step, which can be seen as one step
Fourier transform in this case), λ0 = λ/µ is a ratio of of a gradient descent in the sense of objective minimisation,
regularization parameters from Eq. 4, Ω is an index set of the or a more general approximation function regressing the
(i+1)
− x(i) , i.e. the distance required to move
acquired k-space samples and Λ is a diagonal matrix. Please difference z
refer to [18] for more details of formulating Eq. 4 as a data to the next state. Moreover, note that in every iteration,
consistency layer in a neural network. Eq. 3a is the proximal CRNN updates its internal state H given an input which is
operator of the prior R, and instead of explicitly determining discussed shortly. As such, CRNN also allows information to
the form of the regularisation term, we propose to directly be propagated efficiently across iterations, in contrast to the
learn the proximal operator by using a convolutional recurrent sequential models using CNNs which collapse the intermediate
feature representation to z(i) .
neural network (CRNN).
In order to exploit the dynamic nature and the temporal
Previous deep learning approaches such as Deep-ADMM net
[17] and method proposed by Schlemper et al. [18] unroll the redundancy of our data, we further propose to jointly model the
traditional optimisation algorithm. Hence, their models learn a recurrence evolving over time for dynamic MRI reconstruction.
sequence of transition x(0) → z(1) → x(1) → · · · → z(N ) → The proposed CRNN-MRI network and CRNN block are shown
x(N ) to reconstruct the image, where each state transition at in Fig. 2(a), in which CRNN block comprised of 5 components:
stage (i) is an operation such as convolutions independently
1) bidirectional convolutional recurrent units evolving over
parameterised by θ, nonlinearities or a data consistency step.
time and iterations (BCRNN-t-i),
However, since the network implicitly learns some form of
2) convolutional recurrent units evolving over iterations only
proximal operator at each iteration, it may be redundant to
(CRNN-i),

4

(a) Traditional Optimisation

(c) CRNN-MRI
(𝑖)

argmin

(𝑖)
𝑥𝑟𝑒𝑐

CRNN

𝑥𝑟𝑒𝑐

𝐳

⨁
𝐳

DC

𝑖

DC

𝑖

𝐳

(b) CNN
(0)

𝑥𝑟𝑒𝑐

CNN(1)

𝐳

1

DC

(𝑖−1)

…

CNN(i)

𝑥𝑟𝑒𝑐

𝐳

𝑖

DC

Fig. 1: (a) Traditional optimisation algorithm using variable splitting and alternate minimisation approach, (b) the optimisation
unrolled into a deep convolutional network incorporating the data consistency step, and (c) the proposed architecture which
models optimisation recurrence.
1

𝑖

𝐱 𝑟𝑒𝑐

(a) CRNN-MRI

Ω

DC

⨁

k-space samples

𝐱 𝑟𝑒𝑐

2

𝐱 𝑟𝑒𝑐

DC

DC

DC

⨁

⨁

⨁

(1)
𝐇5

(2)
𝐇5

𝐇5

Unfolding
(iteration)

(1)

𝐇4

(1)

𝐇3

(1)

𝐇2

(1)

𝐇𝟏

𝐇4

CRNN

𝑁

𝐱 𝑟𝑒𝑐

(b)

(2)

(N)

Layer 5
(CNN)

(N)

Layer 4
(CRNN-i)

(N)

Layer 3
(CRNN-i)

(N)

Layer 2
(CRNN-i)

(N)

Layer 1
(BCRNN-t-i)

𝐇4

..

CNN

𝐇3

(2)

𝐇3

..

CRNN-i

𝐇2

(2)

𝐇2

..

CRNN-i
𝐇1

(2)

𝐇1

CRNN-i

𝐱 𝑟𝑒𝑐

𝐇𝑙,𝑡−1)

𝐇𝑙,𝑡

𝑖

(c)

𝑖−1

𝐱 𝑟𝑒𝑐

Feed-forward
convolution

Recurrent
convolution
(iteration)

Recurrent
convolution
(time)

1

𝐱𝑢

BCRNN-t-i

𝑖

𝑖

…

𝐇 𝑙,𝑡−1

…

𝐇𝑙,𝑡−1

𝑖−1

𝑖−1

𝐇𝑙,𝑡

𝑖

𝐇𝑙,𝑡+1

𝑖

𝑖−1

…

𝐇 𝑙,𝑡+1

…

𝑖

𝐇𝑙,𝑡

𝐇𝑙,𝑡+1

𝐇𝑙

𝑖

𝑖

𝐇 𝑙,𝑡

𝑖

𝐇𝑙,𝑡−1

𝑁−1

𝐱 𝑟𝑒𝑐

𝐇𝑙,𝑡+1

𝑖

𝐇𝑙−1,
𝑡−1

𝑖

𝐇𝑙−1,𝑡

𝑖

𝐇𝑙
𝐇𝑙

𝑖

𝑖

𝑖

𝐇𝑙−1,
𝑡+1

Fig. 2: (a) The overall architecture of proposed CRNN-MRI network for MRI reconstruction. (b) The structure of the proposed
(0)
network when unfolded over iterations, in which xrec = xu . (c) The structure of BCRNN-t-i layer when unfolded over the time
sequence. The black arrows indicate feed-forward convolutions. The blue arrows and red arrows indicate recurrent convolutions
over iterations and the time sequence respectively.

3) 2D convolutional neural network (CNN),
4) residual connection and
5) DC layers.
We introduce details of the components of our network in the
following subsections.
1) CRNN-i: As aforementioned, we encapsulate the iterative
optimisation procedures explicitly with RNNs. In the CRNNi unit, the iteration step is viewed as the sequential step in
the vanilla RNN. If the network is unfolded over the iteration
dimension, the network can be illustrated as in Fig. 2(b), where
information is propagated between iterations. Here we use H

to denote the feature representation of our sequence of frames
(i)
throughout the network. Hl denotes the representation at layer
l (subscript) and iteration step i (superscript). Therefore, at
(i−1)
iteration (i), given the input Hl−1 and the previous iteration’s
(i−1)
(i)
hidden state Hl
, the hidden state Hl at layer l of a CRNNi unit can be formulated as:
(i)

(i)

(i−1)

Hl = σ(Wl ∗ Hl−1 + Wi ∗ Hl

+ Bl ).

(7)

Here ∗ represents convolution operation, Wl and Wi
represent the filters of input-to-hidden convolutions and hiddento-hidden recurrent convolutions evolving over iterations

5

(i)

respectively, and Bl represents a bias term. Here Hl is backward direction. Therefore, for the formulation of BCRNNthe representation of the whole T sequence with shape t-i unit, given (1) the current input representation of the l-th
(batchsize, T, nc , Dx , Dy ), where nc is the number of channels layer at time frame t and iteration step i, which is the output
which is 2 at the input and output but is greater while processing representation from (l − 1)-th layer H(i)
l−1,t , (2) the previous
inside the network, and the convolutions are computed on the iteration’s hidden representation within the same layer H(i−1) ,
l,t
→
− (i)
last two dimensions. The latent features are activated by the
(3) the hidden representation of the past time frame H l,t−1 , and
rectifier linear unit (ReLU) as a choice of nonlinearity, i.e.
←
−(i)
the hidden representation of the future time frame H l,t+1 , then
σ(x) = max(0, x).
of the current l-th layer of time
The CRNN-i unit offers several advantages compared to the hidden state representation
(i)
frame
t
at
iteration
i,
H
with
shape (batchsize, nc , Dx , Dy ),
l,t
independently unrolling convolutional filters at each stage.
can
be
formulated
as:
Firstly, compared to CNNs where the latent representation
−
→(i)
−
→
−
→
(i)
(i−1)
from the previous state is not propagated, the hidden-to-hidden
H l,t = σ(Wl ∗ Hl−1,t + Wt ∗ H il,t−1 + Wi ∗ Hl,t + B l ),
←
−(i)
←
−(i)
←
−
iteration connections in CRNN-i units allow contextual spatial
(i)
(i−1)
H l,t = σ(Wl ∗ Hl−1,t + Wt ∗ H l,t+1 + Wi ∗ Hl,t + B l ),
information gathered at previous iterations to be passed to the
−
→(i) ←
−(i)
(i)
Hl,t = H l,t + H l,t ,
future iterations. This enables the reconstruction step at each
(8)
iteration to be optimised not only based on the output image
Similar
to
the
notation
in
Section
III-B1,
W
represents
the
t
but also based on the hidden features from previous iterations,
evolving over time. When
where the hidden connection convolutions can "memorise" the filters of recurrent convolutions
(1)
l
=
1
and
i
=
1,
H
=
x
,
that is the t-th frame of
u
t
0,t
useful features to avoid redundant computation. Secondly, as
undersampled
input
data,
and
when
l = 1 and i = 2, ...T ,
the iteration number increases, the effective receptive field of a
(i)
(i−1)
H
=
x
,
which
stands
for
the
t-th frame of the
rec
t
0,t
CRNN-i unit in the spatial domain also expands whereas CNN
intermediate
reconstruction
result
from
iteration
i − 1. For
resets it at each iteration. This property allows the network to
− (i)
←
−(i)
(0) →
H
,
and
,
they
are
set
to
be
zero
initial
hidden
H
H
l,t
l,0
l,T +1
further improve the reconstruction by allowing it to have better
states.
contextual support. In addition, since the weight parameters
The temporal connections of BCRNN-t-i allow information
are shared across iterations, it greatly reduces the number
to
be propagated across the whole T time frames, enabling it to
of parameters compared to CNNs, potentially offering better
learn
the differences and correlations of successive frames. The
generalization properties.
filter
responses of recurrent convolutions evolving over time
In this work, we use a vanilla RNN to model the recurrence
express dynamic changing biases, which focus on modelling
due to its simplicity. Note this can be naturally generalised
the temporal changes across frames, while the filter responses
to other RNN units, such as long short-term memory (LSTM)
of recurrent convolutions over iterations focus on learning the
and gated recurrent unit (GRU), which are considered to have
spatial refinement across consecutive iteration steps. In addition,
better memory properties, although using these units would
we
note that learning recurrent layers along the temporal
significantly increase computational complexity.
direction is different to using 3D convolution along the space
2) BCRNN-t-i: Dynamic MR images exhibit high temporal and temporal direction. 3D convolution seeks invariant features
redundancy, which is often exploited as a-priori knowledge to across space-time, hence several layers of 3D convolutions are
regularise the reconstruction. Hence, it is also beneficial for the required before the information from the whole sequence can
network to learn the dynamics of sequences. To this extent, we be propagated to a particular time frame. On the other hand,
propose a bidirectional convolutional recurrent unit (BCRNN- learning recurrent 2D convolutions enables the model to easily
t-i) to exploit both temporal and iteration dependencies jointly. and efficiently propagate the information through time, which
BCRNN-t-i includes three convolution layers: one on the input also yields fewer parameters and a lower computational cost.
which comes into the unit from the previous layer, one on
In summary, the set of hidden states for a CRNN block
−(i) →
− (i)
the hidden state from the past and future time frames and
(i)
(i) ←
to update at iteration i is H = {Hl , Hl,t , H l,t , H l,t }, for
the one on the hidden state from the previous iteration of
l = 1, . . . , L and t = 1, . . . , T , where L is the total number of
the unit (Fig. 2(c)). Note that we simultaneously consider
layers in the CRNN block and T is the total number of time
temporal dependencies from past and future time frames, and
frames.
the encoding weights are shared for both directions. The output
for the BCRNN-t-i layer is obtained by summing the feature
maps learned from both directions. The illustration figure of C. Network Learning
the unit when it is unfolded over time sequence is shown in
Given the training data S of input-target pairs (xu , xt ), the
Fig. 2(c).
network learning proceeds by minimizing the pixel-wise mean
As we need to propagate information along temporal squared error (MSE) between the predicted reconstructed MR
dimensions in this unit, here we introduce an additional index t image and the fully sampled ground truth data:
in the notation to represent the variables related with time frame
X
1
2
(i)
L (θ) =
kxt − xrec k2
(9)
t. Here Hl,t represents feature representations at l-th layer,
n
S
→
− (i)
(xu ,xt )∈S
time frame t, and at iteration i, H l,t denotes the representations
calculated when information is propagated forward inside the where θ = {Wl , Wi , Wt , Bl }, l = 1 . . . L, and nS stands
←
−(i)
BCRNN-t-i unit, and similarly, H l,t denotes the one in the for the number of samples in the training set S. Note that

6

the total number of time sequences T and iteration steps N
assumed by the network before performing the reconstruction
is a free parameter that must be specified in advance. The
network weights were initialised using He initialization [26]
and it was trained using the Adam optimiser [27]. During
training, gradients were hard-clipped to the range of [−5, 5]
to mitigate the gradient explosion problem. The network was
implemented in Python using Theano and Lasagne libraries.

aliasing artefact. The evaluation was done via a 3-fold cross
validation. The minibatch size during the training was set to
1, and we observed that the performance can reach a plateau
within 6 × 104 backpropagations.
B. Evaluation Method

We compared the proposed method with the representative
algorithms of the CS-based dynamic MRI reconstruction, such
as k-t FOCUSS [1] and k-t SLR [2], and two variants of
IV. E XPERIMENTS
3D CNN networks named 3D CNN-S and 3D CNN in our
experiments. The built baseline 3D CNN networks share the
A. Dataset and Implementation Details
same architecture with the proposed CRNN-MRI network but
The proposed method was evaluated using a complexall the recurrent units and 2D CNN units were replaced with
valued MR dataset consisting of 10 fully sampled short-axis
3D convolutional units, that is, in each iteration, the 3D CNN
cardiac cine MR scans. Each scan contains a single slice
block contain 5 layers of 3D convolutions, one DC layer and a
SSFP acquisition with 30 temporal frames, which have a
residual connection. Here 3D CNN-S refers to network sharing
320 × 320 mm field of view and 10 mm thickness. The
weights across iterations, however, this does not employ the
raw data consists of 32-channel data with sampling matrix
hidden-to-hidden connection as in the CRNN-i unit. The 3D
size 192 × 190, which was then zero-filled to the matrix size
CNN-S architecture was chosen so as to make a fair comparison
256 × 256. The raw multi-coil data was reconstructed using
with the proposed model using a comparable number of network
SENSE [28] with no undersampling and retrospective gating.
parameters. In contrast, 3D CNN refers to the network without
Coil sensitivity maps were normalized to a body coil image
weight sharing, in which the network capacity is N = 10 times
and used to produce a single complex-valued reconstructed
of that of 3D CNN-S, and approximately 12 times more than
image. In experiments, the complex valued images were backthat of our first proposed method (Proposed-A).
transformed to regenerate k-space samples, simulating a fully
Reconstruction results were evaluated based on the followsampled single-coil acquisition. The input undersampled image
ing quantitative metrics: MSE, peak-to-noise-ratio (PSNR),
sequences were generated by randomly undersampling the
structural similarity index (SSIM) [29] and high frequency
k-space samples using Cartesian undersampling masks, with
error norm (HFEN) [30]. The choice of the these metrics was
undersampling patterns adopted from [1]: for each frame the
made to evaluate the reconstruction results with complimentary
eight lowest spatial frequencies were acquired, and the sampling
emphasis. MSE and PSNR were chosen to evaluate the overall
probability of k-space lines along the phase-encoding direction
accuracy of the reconstruction quality. SSIM put emphasis
was determined by a zero-mean Gaussian distribution. Note
on image quality perception. HFEN was used to quantify the
that the undersampling rates are stated with respect to the
quality of the fine features and edges in the reconstructions,
matrix size of raw data, which is 192 × 190.
and here we employed the same filter specification as in [30],
The architecture of the proposed network used in the
[31] with the filter kernel size 15 × 15 pixels and a standard
experiment is shown in Fig. 2: each iteration of the CRNN
deviation of 1.5 pixels. For PSNR and SSIM, it is the higher
block contains five units: one layer of BCRNN-t-i, followed
the better, while for MSE and HFEN, it is the lower the better.
by three layers of CRNN-i units, and followed by a CNN
unit. For all CRNN-i and BCRNN-t-i units, we used a kernel
size k = 3 and the number of filters was set to nf = 64 for C. Results
Proposed-A and nf = 128 for Proposed-B in Table I. The CNN
The comparison results of all methods are reported in Table I,
after the CRNN-i units contains one convolution layer with where we evaluated the quantitative metrics, network capacity
k = 3 and nf = 2, which projects the extracted representation and reconstruction time. Numbers shown in Table I are mean
back to the image domain which contains complex-valued values of corresponding metrics with standard deviation of
images expressed using two channels. The output of the CRNN different subjects in parenthesis. Bold numbers in Table I
block is connected to the residual connection, which sums indicate the better performance of the proposed methods than
the output of the block with its input. Finally, we used DC the competing ones. Compared with the baseline method (k-t
layers on top of the CRNN output layers. During training, the FOCUSS and k-t SLR), the proposed methods outperform them
iteration step is set to be N = 10, and the time sequence by a considerable margin at different acceleration rates. When
for training is T = 30. Note that this architecture is by no compared with deep learning methods, note that the network
means optimal and more layers can be added to increase the capacity of Proposed-A is comparable with that of 3D CNN-S
ability of our network to better capture the data structures. and the capacity of Propose-B is around one third of that of 3D
While the original sequence has size 256 × 256 × T , we extract CNN. Though their capacities are much smaller, both Proposedpatches of size 256 × Dpatch × T , where Dpatch is the patch A and Proposed-B outperform 3D CNN-S and 3D CNN for
size and the direction of patch extraction corresponds to the all acceleration rates by a large margin, which shows the
frequency-encoding direction. Note that since we only consider competitiveness and effectiveness of our method. In addition,
Cartesian undersampling, the aliasing occurs only along the we can see a substantial improvement of the reconstruction
phase encoding direction, so patch extraction does not alter the results on all acceleration rates and in all metrics when the

7

TABLE I: Performance comparisons (MSE, PSNR:dB, SSIM, and HFEN) on dynamic cardiac data with different acceleration
rates. MSE is scaled to 10−3 . The bold numbers are better results of the proposed methods than that of the other methods.
Method

k-t FOCUSS

k-t SLR

3D CNN-S

3D CNN

Proposed-A

Proposed-B

-

-

338,946

3,389,460

262,020

1,040,132

6×

MSE
PSNR
SSIM
HFEN

0.592 (0.199)
32.506 (1.516)
0.953 (0.040)
0.211 (0.021)

0.371( 0.155)
34.632 (1.761)
0.970 (0.033)
0.161 (0.016)

0.385 (0.124)
34.370 (1.526)
0.976 (0.008)
0.170 (0.009)

0.275 (0.096)
35.841 (1.470)
0.983 (0.005)
0.138 (0.013)

0.261 (0.097)
36.096 (1.539)
0.985 (0.004)
0.131 (0.013)

0.201 (0.074)
37.230 (1.559)
0.988 (0.003)
0.112 (0.010)

9×

MSE
PSNR
SSIM
HFEN

1.234 (0.801)
29.721 (2.339)
0.922 (0.043)
0.310(0.041)

0.846 (0.572)
31.409 (2.404)
0.951 (0.025)
0.260 (0.034)

0.929 (0.474)
30.838 (2.246)
0.950 (0.016)
0.280 (0.034)

0.605 (0.324)
32.694 (2.179)
0.968 (0.010)
0.215 (0.021)

0.516 (0.255)
33.281 (1.912)
0.972 (0.009)
0.201 (0.025)

0.405 (0.206)
33.379 (2.017)
0.979 (0.007)
0.173 (0.021)

11×

MSE
PSNR
SSIM
HFEN

1.909 (0.828)
27.593 (2.038)
0.880 (0.060)
0.390 (0.023)

1.237 (0.620)
29.577 (2.211)
0.924 (0.034)
0.327 (0.028)

1.472 (0.733)
28.803 (2.151)
0.925 (0.022)
0.363 (0.041)

0.742 (0.325)
31.695 (1.985)
0.960 (0.010)
0.257 (0.029)

0.688 (0.290)
31.986 (1.885)
0.964 (0.009)
0.248 (0.033)

0.610 (0.300)
32.575 (1.987)
0.968 (0.011)
0.227 (0.030)

15s

451s

8s

8s

3s

6s

Capacity

Time

images and their corresponding error maps from different
reconstruction methods. As one can see, our proposed model
(Proposed-B) can produce more faithful reconstructions for
those parts of the image around the myocardium where there
are large temporal changes. This is reflected by the fact that
RNNs effectively use a larger receptive field to capture the
characteristics of aliasing seen within the anatomy. For the
3D CNN approaches, it is also observed that it is not able to
denoise the background region. This could be explained by the
fact that 3D CNN only exploits local information due to the
small filter sizes it used, while in contrast, the proposed CRNN
improves the denoising of the background region because of its
larger receptive field sizes. Their temporal profiles at x = 120
are shown in Fig. 5. Similarly, one can see that the proposed
model has overall much smaller error, faithfully modelling
the dynamic data. This suggests its capability to learn motion
Fig. 3: Mean PSNR values (Proposed-B) vary with the number compensation implicitly between frames although the network
of iterations at test time on data with different acceleration is trained for the dynamic image reconstruction. It could be
rates. Here ac stands for acceleration rate.
due to the fact that spatial and temporal features are learned
separately in the proposed model while 3D CNN seeks invariant
feature learning across space and time.
number of network parameters is increased for the proposed
In terms of speed, the proposed RNN-based reconstruction
method (Proposed-B), and therefore we will only show the is faster than the 3D CNN approaches because it only
results from Proposed-B in the following. The number of performs convolution along time once per iteration, removing
iterations used by the network at test time is set to be the the redundant 3D convolutions which are computationally
same as the training stage, which is N = 10, however, if expensive. Reconstruction time of 3D CNN and the proposed
the iteration number is increased up to N = 17, it shows an methods reported in Table I were calculated on a GPU GeForce
improvement of 0.324dB on average. Fig. 3 shows the model’s GTX 1080, and the time for k-t FOCUSS and k-t SLR were
performance varying with the number of iterations at test time. calculated on CPU.
In fact, for accelerated imaging, higher undersampling factors
significantly add aliasing to the initial zero-filled reconstruction,
V. D ISCUSSION AND CONCLUSION
making the reconstruction more challenging. This suggests that
while the 3D CNN possesses higher modelling capacity owing
In this work, we have demonstrated that the presented
to its large number of parameters, it may not necessarily be network is capable of producing faithful image reconstructions
an ideal architecture to perform dynamic MR reconstruction, from highly undersampled data, both in terms of various
presumably because the simple CNN is not as efficient as quantitative metrics as well as inspection of error maps.
propagating the information across the whole sequence.
In contrast to unrolled deep network architectures proposed
A comparison of the visualization results of a reconstruction previously, we modelled the recurrent nature of the optimisation
from 9× acceleration is shown in Fig. 4 with the reconstructed iteration using hidden representations with the ability to

8

Fig. 4: The comparison of reconstructions on spatial dimension with their error maps. (a) Ground Truth (b) Undersampled
image by acceleration factor 9 (c,d) Proposed-B (e,f) 3D CNN (g,h) 3D CNN-S (i,j) k-t FOCUSS (k,l) k-t SLR

Additionally, current analysis only considers a single coil setup.
In the future, we will also aim at investigating such methods in
a scenario where multiple coil data from parallel MR imaging
can be used jointly for higher acceleration acquisition.
To conclude, inspired by variable splitting and alternate
minimisation strategies, we have presented an end-to-end deep
learning solution, CRNN-MRI, for accelerated dynamic MRI
reconstruction, with a forward, CRNN block implicitly learning
iterative denoising interleaved by data consistency layers to
Fig. 5: The comparison of reconstructions along temporal enforce data fidelity. In particular, the CRNN architecture is
dimension with their error maps. (a) Ground Truth (b) Under- composed of the proposed novel variants of convolutional
sampled image by acceleration factor 9 (c,d) Proposed-B (e,f) recurrent unit which evolves over two dimensions: time and
iterations. The proposed network is able to learn both the
3D CNN (g,h) 3D CNN-S (i,j) k-t FOCUSS (k,l) k-t SLR
temporal dependency and the iterative reconstruction process
effectively, and outperformed the other competing methods in
retain and propagate information across the optimisation steps. terms of computational complexity, reconstruction accuracy
Compared with 3D CNN models, the proposed methods have a and speed for different undersampling rates.
much lower network capacity but still have a higher accuracy,
reflecting the effectiveness of our architecture. This is due to
R EFERENCES
the ability of the proposed RNN units to increase the receptive
[1] H. Jung, J. C. Ye, and E. Y. Kim, “Improved k–t BLAST and k–t SENSE
field size while iteration steps increase, as well as to efficiently
using FOCUSS,” Physics in medicine and biology, vol. 52, no. 11, p.
propagate information across the temporal direction. In addition,
3201, 2007.
our network also offers very fast reconstruction on a GPU
[2] S. G. Lingala, Y. Hu, E. Dibella, and M. Jacob, “Accelerated dynamic
MRI exploiting sparsity and low-rank structure: K-t SLR,” IEEE
GeForce GTX 1080 compared with other competing methods.
Transactions on Medical Imaging, vol. 30, no. 5, pp. 1042–1054, 2011.
In this work, we modeled the recurrence using the relatively
[3] J. Schlemper, J. Caballero, J. V. Hajnal, A. Price, and D. Rueckert, “A
simple (vanilla) RNN architecture. For the future work, we will
deep cascade of convolutional neural networks for dynamic MR image
reconstruction,” arXiv preprint arXiv:1704.02422, 2017.
explore other recurrent units such as LSTM or GRU. As they are
[4] D. Batenkov, Y. Romano, and M. Elad, “On the global-local dichotomy
trained to explicitly select what to remember, they may allow
in sparsity modeling,” arXiv preprint arXiv:1702.03446, 2017.
the units to better control the flow of information and could
[5] J. Tsao, P. Boesiger, and K. P. Pruessmann, “k-t BLAST and k-t
SENSE: Dynamic MRI with high frame rate exploiting spatiotemporal
reduce the number of iterations required for the network to
correlations,” Magnetic Resonance in Medicine, vol. 50, no. 5, pp. 1031–
generate high-quality output. In addition, we have found that the
1042, 2003.
majority of errors between the reconstructed image and the fully
[6] J. Caballero, A. N. Price, D. Rueckert, and J. V. Hajnal, “Dictionary
learning and time sparsity for dynamic MR data reconstruction,” IEEE
sampled image lie at the part where motion exists, indicating
Transactions on Medical Imaging, vol. 33, no. 4, pp. 979–994, 2014.
that motion exhibits a challenge for such dynamic sequence
[7] R. Otazo, E. Candès, and D. K. Sodickson, “Low-rank plus sparse
reconstruction, and CNNs or RNNs trained for reconstruction
matrix decomposition for accelerated dynamic MRI with separation of
background and dynamic components,” Magnetic Resonance in Medicine,
loss cannot fully learn such motion compensation implicitly
vol. 73, no. 3, pp. 1125–1136, 2015.
during training. Thus it will be interesting to explore the
[8] M. Lustig, D. L. Donoho, J. M. Santos, and J. M. Pauly, “Compressed
benefits of doing simultaneous motion compensation and image
sensing mri,” IEEE signal processing magazine, vol. 25, no. 2, pp. 72–82,
reconstruction for the improvement in the dynamic region.
2008.

9

[9] K. Hammernik, T. Klatzer, E. Kobler, M. P. Recht, D. K. Sodickson,
T. Pock, and F. Knoll, “Learning a variational network for reconstruction
of accelerated MRI data,” arXiv preprint arXiv:1704.00447, 2017.
[10] C. Dong, C. C. Loy, K. He, and X. Tang, “Learning a deep convolutional network for image super-resolution,” in European Conference on
Computer Vision. Springer, 2014, pp. 184–199.
[11] O. Ronneberger, P. Fischer, and T. Brox, “U-Net: Convolutional Networks
for Biomedical Image Segmentation,” Medical Image Computing and
Computer-Assisted Intervention, pp. 234–241, 2015.
[12] D. Lee, J. Yoo, and J. C. Ye, “Deep artifact learning for compressed
sensing and parallel MRI,” arXiv preprint arXiv:1703.01120, 2017.
[13] Y. S. Han, J. Yoo, and J. C. Ye, “Deep learning with domain adaptation for accelerated projection reconstruction MR,” arXiv preprint
arXiv:1703.01135, 2017.
[14] S. Wang, Z. Su, L. Ying, X. Peng, and D. Liang, “Exploiting deep
convolutional neural network for fast magnetic resonance imaging,” in
ISMRM 24th Annual Meeting and Exhibition, 2016.
[15] S. Wang, N. Huang, T. Zhao, Y. Yang, L. Ying, and D. Liang, “1D Partial
Fourier Parallel MR imaging with deep convolutional neural network,”
in ISMRM 25th Annual Meeting and Exhibition, vol. 47, no. 6, 2017,
pp. 2016–2017.
[16] J. Adler and O. Öktem, “Learned primal-dual reconstruction,” arXiv
preprint arXiv:1707.06474, 2017.
[17] J. Sun, H. Li, Z. Xu et al., “Deep ADMM-Net for compressive sensing
mri,” in Advances in Neural Information Processing Systems, 2016, pp.
10–18.
[18] J. Schlemper, J. Caballero, J. V. Hajnal, A. Price, and D. Rueckert,
“A deep cascade of convolutional neural networks for MR image
reconstruction,” arXiv preprint arXiv:1703.00555, 2017.
[19] J. Adler and O. Öktem, “Solving ill-posed inverse problems using iterative
deep neural networks,” arXiv preprint arXiv:1704.04058, 2017.
[20] K. Zhang, W. Zuo, S. Gu, and L. Zhang, “Learning deep cnn denoiser
prior for image restoration,” arXiv preprint arXiv:1704.03264, 2017.
[21] J. Chang, C.-L. Li, B. Póczos, B. Kumar, and A. C. Sankaranarayanan,
“One network to solve them all—solving linear inverse problems using
deep projection models,” arXiv preprint arXiv:1703.09912, 2017.
[22] K. Gregor, I. Danihelka, A. Graves, D. Rezende, and D. Wierstra, “Draw:
A recurrent neural network for image generation,” in Proceedings of
The 32nd International Conference on Machine Learning, 2015, pp.
1462–1471.
[23] M. Liang and X. Hu, “Recurrent convolutional neural network for object
recognition,” in Proceedings of the IEEE Conference on Computer Vision
and Pattern Recognition, 2015, pp. 3367–3375.
[24] J. Kuen, Z. Wang, and G. Wang, “Recurrent attentional networks for
saliency detection,” in Proceedings of the IEEE Conference on Computer
Vision and Pattern Recognition, 2016, pp. 3668–3677.
[25] Y. Huang, W. Wang, and L. Wang, “Bidirectional recurrent convolutional
networks for multi-frame super-resolution,” in Advances in Neural
Information Processing Systems, 2015, pp. 235–243.
[26] K. He, X. Zhang, S. Ren, and J. Sun, “Delving deep into rectifiers:
Surpassing human-level performance on imagenet classification,” in
Proceedings of the IEEE international conference on computer vision,
2015, pp. 1026–1034.
[27] D. Kingma and J. Ba, “Adam: A method for stochastic optimization,”
arXiv preprint arXiv:1412.6980, 2014.
[28] K. P. Pruessmann, M. Weiger, M. B. Scheidegger, P. Boesiger et al.,
“SENSE: sensitivity encoding for fast MRI,” Magnetic resonance in
medicine, vol. 42, no. 5, pp. 952–962, 1999.
[29] Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, “Image
quality assessment: from error visibility to structural similarity,” IEEE
transactions on image processing, vol. 13, no. 4, pp. 600–612, 2004.
[30] S. Ravishankar and Y. Bresler, “Mr image reconstruction from highly
undersampled k-space data by dictionary learning,” IEEE transactions
on medical imaging, vol. 30, no. 5, pp. 1028–1041, 2011.
[31] X. Miao, S. G. Lingala, Y. Guo, T. Jao, M. Usman, C. Prieto, and K. S.
Nayak, “Accelerated cardiac cine mri using locally low rank and finite
difference constraints,” Magnetic resonance imaging, vol. 34, no. 6, pp.
707–714, 2016.

