{'fraction_non_alphanumeric': 0.04013377926421405, 'fraction_numerical': 0.011705685618729096, 'mean_word_length': 4.667409967666407, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Online Adaptive Image Reconstruction (OnAIR) Using Dictionary Models

Sparsity and low-rank models have been popular for reconstructing images and videos from limited or corrupted measurements. Dictionary or transform learning methods are useful in applications such as denoising, inpainting, and medical image reconstruction. This paper proposes a framework for online (or time-sequential) adaptive reconstruction of dynamic image sequences from linear (typically undersampled) measurements. We model the spatiotemporal patches of the underlying dynamic image sequence as sparse in a dictionary, and we simultaneously estimate the dictionary and the images sequentially from streaming measurements. Multiple constraints on the adapted dictionary are also considered such as a unitary matrix, or low-rank dictionary atoms that provide additional efficiency or robustness. The proposed online algorithms are memory efficient and involve simple updates of the dictionary atoms, sparse coefficients, and images. Numerical experiments demonstrate the usefulness of the proposed methods in inverse problems such as video reconstruction or inpainting from noisy, subsampled pixels, and dynamic magnetic resonance image reconstruction from very limited measurements.


I. INTRODUCTION
Models of signals and images based on sparsity, low-rank, and other properties are useful in image and video processing. In ill-posed or ill-conditioned inverse problems, it is often useful to employ signal models that reflect known or assumed properties of the latent images. Such models are often used to construct appropriate regularization. For example, the sparsity of images in wavelet or discrete cosine transform (DCT) domains has been exploited for image and video reconstruction tasks [1]- [3]. In particular, the learning of such models has been explored in various settings [4]- [10], where they may potentially outperform fixed models since they adapt to signals and signal classes.
There has been growing interest in such dictionary learningbased image restoration or reconstruction methods [11]- [13]. For example, in blind compressed sensing [10], [12], [13], a dictionary for the underlying image or video is estimated together with the image from undersampled measurements. This allows the dictionary to adapt to the current data, which may enable learning novel features, providing improved reconstructions.
For inverse problems involving large-scale or streaming data, e.g., in interventional imaging or restoring (e.g, denoising, inpainting) large or streaming videos, etc., it is often critical to obtain reconstructions in an online or time-sequential (or data-sequential) manner to limit latency. Batch methods that process all the data at once are often prohibitively expensive in terms of time and memory usage, and infeasible for streaming data. Methods for online learning of dictionary and sparsifying transform models from streaming signals including noisy data have been recently proposed [14]- [18] and shown to outperform state-of-the-art methods [19]- [21] including deep learning methods for video denoising.
This paper focuses on methods for dynamic image reconstruction from limited measurements such as video or medical imaging data. One such important class of dynamic image reconstruction problems arises in dynamic magnetic resonance imaging (dMRI), where the data are inherently or naturally undersampled because the object is changing as the data (samples in k-space or Fourier space of the object acquired sequentially over time) is collected. Various techniques have been proposed for reconstructing such dynamic image sequences from limited measurements [22]- [24]. Such methods may achieve improved spatial or temporal resolution by using more explicit signal models compared to conventional approaches (such as k-space data sharing in dMRI, where data is pooled in time to make sets of k-space data with sufficient samples, such as in the form of a Casorati matrix [25]); these methods typically achieve increased accuracy at the price of increased computation.
While some reconstruction techniques are driven by sparsity models and assume that the image sequence is sparse in some transform domain or dictionary [22], other methods exploit low-rank or other kinds of sophisticated models [24]- [30]. For example, L+S methods [24], [30] assume that the image sequence can be decomposed as the sum of low-rank and sparse (either directly sparse or sparse in some known transform domain) components that are estimated from measurements. Dictionary learning-based approaches including with low-rank models provide promising performance for dynamic image reconstruction [12], [31], [32]. Although these methods allow adaptivity to data and provide improved reconstructions, they involve batch processing and typically expensive computation and memory use. Next, we outline our contributions, particularly a new online framework that address these issues of state-of-the-art batch reconstruction methods.

A. Contributions
This paper investigates a framework for Online Adaptive Image Reconstruction, dubbed OnAIR, exploiting learned arXiv:1809.01817v1 [stat.ML] 6 Sep 2018 dictionaries. In particular, we model spatiotemporal patches of the underlying dynamic image sequence as sparse in an (unknown) dictionary, and propose a method to jointly and sequentially (over time) estimate the dictionary, sparse codes, and images from streaming measurements. Various constraints are considered for the dictionary model such as a unitary matrix, or a matrix whose atoms/columns are low-rank when reshaped into spatio-temporal matrices. The proposed OnAIR algorithms involve simple and efficient updates and require a small, fixed amount of data to be stored in memory at a time, greatly reducing the computation and memory demands during processing compared to conventional batch methods that process all of the data simultaneously.
Numerical experiments demonstrate the effectiveness of the proposed OnAIR methods for performing video inpainting from subsampled and noisy pixels, and dynamic MRI reconstruction from very limited k-t space measurements. The experiments show that our proposed methods are able to learn dictionaries adaptively from corrupted measurements with important representational features that improve the quality of the reconstructions compared to non-adaptive schemes. Moreover, the OnAIR methods provide better image quality than batch schemes with a single adapted dictionary since they allow the dictionary to change over time for changing (dynamic) videos. At the same time, they require much lower runtimes than adaptive batch (or offline) reconstructions.
Short versions of this work appeared recently in [33] and [34]. This paper extends those initial works by exploring OnAIR methods with multiple constraints on the learned models including a novel Unitary Dictionary (OnAIR-UD) constraint, and a Low-rank Dictionary atom (OnAIR-LD) constraint. Such constraints may offer a degree of algorithmic efficiency and robustness to artifacts in limited data settings. For example, the updates per time instant in OnAIR-UD are simpler and non-sequential, enabling it to be faster in practice than the other OnAIR schemes. Importantly, this paper reports extensive numerical experiments investigating and evaluating the performance of the proposed online methods in multiple inverse problems, and comparing their performance against recent related online and batch methods in the literature. Finally, we also compare the performance of the OnAIR methods to an oracle online scheme, where the dictionary is learned offline from the ground truth data, and show that the proposed online learning from highly limited and corrupted streaming measurements performs comparably to the oracle scheme.

B. Organization
The rest of this paper is organized as follows. Section II reviews the dictionary learning framework that forms the basis of our online schemes. Section III describes our formulations for online adaptive image reconstruction and Section IV presents algorithms for solving the problems. Section V presents extensive numerical experiments that demonstrate the promising performance of our proposed methods in inverse problem settings such as video inpainting and dynamic MRI reconstruction from highly limited data. Finally, Section VI concludes with proposals for future work.

II. DICTIONARY LEARNING MODELS
Here, we briefly review some prior works on dictionary learning that helps build our OnAIR framework. Given a set of signals (or vectorized image patches) that are represented as columns of a matrix P ∈ C n×M , the goal of dictionary learning (DL) is to learn a dictionary D ∈ C n×m and a matrix Z ∈ C m×M of sparse codes such that P ≈ DZ. Traditionally, the DL problem is often formulated [4] as follows: where z i and d j denote the ith column of Z and the jth column of D, respectively, and s denotes a target sparsity level for each signal. Here, the 0 "norm" counts the number of non-zero entries of its argument, and the columns of D are set to unit norm to avoid scaling ambiguity between D and Z [35]. Alternatives to (1) exist that replace the 0 "norm" constraint with other sparsity-promoting constraints, or enforce additional properties on the dictionary [36]- [38], or enable online dictionary learning [14].
Dictionary learning algorithms [4], [14], [39]- [42] typically attempt to solve (1) or its variants in an alternating manner by performing a sparse coding step (updating Z) followed by a dictionary update step (updating D). Methods such as K-SVD [4] also partially update the coefficients in Z in the dictionary update step, while some recent methods update the variables jointly and iteratively [43]. Algorithms for (1) typically repeatedly update Z (involves an NP-hard sparse coding problem) and tend to be computationally expensive.
The DINO-KAT learning problem [13] is an alternative dictionary learning framework that imposes a low-rank constraint on reshaped dictionary atoms (columns). The corresponding problem formulation is where the 0 "norm" counts the total number of nonzeros in Z. The operator R(·) reshapes dictionary atoms d i ∈ C n into matrices of size n 1 × n 2 for some n 1 and n 2 such that n = n 1 n 2 , and r > 0 is the maximum allowed rank for each reshaped atom. The dimensions of the reshaped atoms can be chosen on an application-specific basis. For example, when learning D from 2D image patches, the reshaped atoms can have the dimensions of the 2D patch. In the case where spatiotemporal (3D) patches are vectorized and extracted from dynamic data, the atoms could be reshaped into space-time (2D) matrices. Spatiotemporal patches of videos often have high temporal correlation, so they may be well represented by a dictionary with low-rank space-time (reshaped) atoms [13]. The parameter λ > 0 in (2) controls the overall sparsity of Z. Penalizing the overall or aggregate sparsity of Z enables variable sparsity levels across training signals (a more flexible model than (1)). The ∞ constraints for L > 0 prevent pathologies (e.g., unbounded algorithm iterates) due to the non-coercive objective [44]. In practice, we set L very large, and the constraint is typically inactive. We proposed an efficient algorithm in [13] for Problem (2) that requires less computation than algorithms such as K-SVD. Another alternative to the DL problems in (1) and (2) involves replacing the constraints in (2) with the unitary constraint D H D = I, where I is the n × n identity matrix. Learned unitary operators work well in image denoising and more general reconstruction problems [45], [46]. Moreover, algorithms for learning unitary dictionaries tend to be computationally cheap [45]. Our OnAIR framework exploits some of the aforementioned dictionary structures.

III. PROBLEM FORMULATIONS
This section formulates data-driven online image reconstruction. First, we propose an online image reconstruction framework based on an adaptive dictionary regularizer as in (2). Let {f t ∈ C Nx×Ny } denote the sequence of dynamic image frames to be reconstructed. We assume that noisy, undersampled linear measurements of these frames are observed. We process the streaming measurements in minibatches, with each minibatch containing measurements ofM ≥ 1 consecutive frames. Let x t denote the vectorized version of the 3D tensor obtained by (temporally) concatenatingM consecutive frames of the unknown dynamic images. In practice, we construct the sequence {x t } using a sliding window (over time) strategy, which may involve overlapping or non-overlapping minibatches. We model the spatiotemporal (3D) patches of each x t as sparse with respect to a latent dictionary D. Under this model, we propose to solve the following online dictionary learning-driven image reconstruction problem for each time t = 1, 2, 3, . . .
Here j indexes time, and y t denotes the (typically undersampled) measurements that are related to the underlying frames x t (that we would like to reconstruct) through the linear sensing operator A t . For example, in video inpainting, A t samples a subset of pixels in x t , or in dynamic MRI, A t corresponds to an undersampled Fourier encoding. The operator P l is a patch extraction matrix that extracts the lth n x × n y × n t spatiotemporal patch from x t as a vector. A total of M (possibly) overlapping 3D patches are assumed. Matrix D ∈ C n×m with n = n x n y n t is the synthesis dictionary to be learned and z t l ∈ C m is the unknown sparse code for the lth patch of x t , with P l x t ≈ Dz t l . Matrix Z t has z t l as its columns, and its aggregate (total) sparsity is penalized. The weights λ S , λ Z ≥ 0 are regularization parameters that control the relative adaptive dictionary regularization and sparsity of Z, respectively, in the model. Problem (P1) has a form often used for online optimization [14], [47]. The objective is a weighted average of the function h(y t , D, with the constraints defining the feasible sets for D and the sparse codes. Problem (P1) jointly estimates the adaptive dictionary model for the patches of x t together with the underlying image frames. Note that, for each time index t, we solve (P1) only for the latest group of frames x t and the latest sparse coefficients Z t , while the previous images and sparse coefficients are set to their estimates from previous minibatches (i.e., x j =x j and Z j =Ẑ j for j < t). However, the dictionary D is adapted to all the spatiotemporal patches observed up to time t. 1 An exponential forgetting factor ρ t−j with 0 < ρ < 1 is used for the terms in (P1), and K t = t j=1 ρ t−j is a normalization constant for the objective. The forgetting factor ρ diminishes the influence of "old" data on the dictionary adaptation process. When the dynamic object or scene changes slowly over time, a large ρ (close to 1) is preferable so that past information has more influence and vice versa.
As written in (P1), the dictionary D is updated based on patches from all previous times; however, the proposed algorithm does not need to store this information during optimization. Indeed, our algorithm in Section IV computes only a few constant-sized matrices that contain the necessary cumulative (over time) information to solve (P1).
When minibatches x t and x t+1 do not overlap (i.e., no common frames), each frame f t is reconstructed exactly once in its corresponding window in (P1). However, it is often beneficial to construct the x t 's using an overlapping sliding window strategy [17], in which case a frame f t may be reconstructed in multiple windows (minibatches of frames). In this case, we independently produce estimatesx t for each time index as indicated in (P1), and then we produce a final estimate of the underlying frame f t by computing a weighted average of the reconstructions of that frame from each window in which it appeared. We found empirically that an exponentially ρ-weighted average (similar to that in (P1)) performed better than alternatives such as an unweighted average of the estimates from each window or using only the most recent reconstruction from the latest window.
In (P1), we imposed a low-rank constraint on the dictionary atoms. As an alternative, we consider constraining the dictionary to be a unitary matrix. The online optimization problem in this case is as follows: where all terms in (P2) are defined as in (P1). Note that (P2) does not require the ∞ -norm constraints on the sparse coefficients Z t because the unitary constraint on the dictionary precludes the possibility of repeated dictionary atoms, which was the motivation for including these constraints in (P1) [44].
IV. ALGORITHMS AND PROPERTIES This section presents the algorithms for Problems (P1) and (P2) and their properties. We propose an alternating minimization-type scheme for (P1) and (P2) and exploit the online nature of the optimization to minimize the costs efficiently. At each time index t, we alternate a few times between updating (D, Z t ) while holding x t fixed (the dictionary learning step) and then updating x t with (D, Z t ) held fixed (the image update step). For each t, we use a warm start for (initializing) the alternating approach. We initialize the dictionary D with the most recent dictionary (D t−1 ).
Frames of x t that were estimated in the previous (temporal) windows are initialized with the most recent ρ-weighted reconstructions, and new frames are initialized using simple approaches (e.g., interpolation in the case of inpainting). Initializing the sparse coefficients Z t with the codes estimated in the preceding window (Ẑ t−1 ) worked well in practice. All updates are performed efficiently and with modest memory usage as will be shown next. Figure 1 provides a graphical flowchart depicting our proposed online scheme. The next subsections present the alternating scheme for each time t in more detail.

A. Dictionary Learning Step for (P1)
Let C t (Z t ) H . Minimizing (P1) with respect to (D, C t ) yields the following optimization problem: is the matrix whose columns contain the patches P l x j for 1 ≤ l ≤ M , and c t i is the ith column of C t . We use a block coordinate descent approach [44] (with few iterations) to update the sparse coefficients c t i and atoms d i (columns of D) in (3) sequentially. For each 1 ≤ i ≤ m, we first minimize (3) with respect to c t i keeping the other variables fixed (sparse coding step), and then we update d i keeping the other variables fixed (dictionary atom update step). These updates are performed in an efficient online manner as detailed next.

1) Sparse Coding
Step: Minimizing (3) with respect to c t i leads to the following subproblem: min where the matrix is defined based on the most recent estimates of the other atoms and sparse coefficients. The solution to (4), assuming L > λ Z , is given by [44]  x t y y t Fig. 1: Flowchart of the proposed online adaptive dictionary learning-driven image reconstruction scheme at time t. The input is a vector y t containing the streaming measurements for the current minibatch of frames,x t denotes the corresponding reconstructed minibatch, and x t denotes intermediate reconstructions in the alternating scheme. In the dictionary learning step (D, C t ) are updated with x t held fixed, by performing block coordinate descent over the columns of C t (sparse coding) and the columns of D (dictionary atom update). Then, the frames x t are updated (the image update step) with (D, C t ) held fixed. This process is repeated a few times, and the final frame estimatesx t are integrated into the streaming reconstruction,x. x t 0 denotes the initialization (warm start).
where H λ Z (·) is the elementwise hard thresholding operator that sets entries with (complex) magnitude less than λ Z to zero and leaves other entries unaffected, 1 M is a length-M vector of ones, and min(·, ·) denote elementwise multiplication and elementwise minimum respectively, and e j∠· is computed elementwise, with ∠ denoting the phase. We do not construct E t i in (6) explicitly; rather we efficiently compute the matrixvector product i based on the most recent estimates of each quantity using sparse matrix-vector operations [44].
2) Dictionary Atom Update Step: Here, we minimize (3) with respect to d i . This update uses past information via the forgetting factor ρ. LetP j := ρ t−j P j andC j := ρ t−j C j denote the ρ-weighted patches and sparse coefficients, respectively, and letP 1:t andC 1:t denote the matrices formed by stacking theP j 's horizontally andC j 's vertically, respectively, for times 1 to t. Finally, defineẼ 1:t i :=P 1:t − k =i d k (c 1:t k ) H using the most recent estimates of all variables, withc 1:t k denoting the kth column ofC 1:t . Using this notation, the minimization of (3) with respect to d i becomes Let U r Σ r V H r be the rank-r truncated singular value decomposition (SVD) of the matrix R(Ẽ 1:t ic 1:t i ) that is obtained by computing the r leading singular vectors and singular values of the full SVD R(Ẽ 1:t Then a solution to (7) is given by [13] R where W is any matrix 2 of appropriate dimension with rank at most r such that W F = 1.
The main computation in (8) is computingẼ 1:t ic 1:t i , since the SVD of the small n y n x × n t (i.e., space-time matrix with n t < n y n x typically) matrix R(Ẽ 1:t ic 1:t i ) has negligible computational cost. In principle, the matrix-vector multiplicatioñ E 1:t ic 1:t i depends on all past information processed by the streaming algorithm; however, it can be recursively computed using constant time and memory. Indeed, observe that where [z] i denotes the ith element of a vector z. The vectors q t i and g t i depend on all previous data, but they can be recursively computed over time as where for each column index i, the matrix C t is understood to contain the latest versions of the sparse codes already updated (sequentially) during the dictionary learning step. Using these recursive formulas, the productẼ 1:t ic 1:t i can be readily computed each time in our algorithm. Thus, the update in (8) can be performed in a fully online manner (i.e., without storing all the past data).
Alternatively, we can collect the vectors q t−1 i and g t−1 i as columns of matrices Q t−1 ∈ C n×m and G t−1 ∈ C m×m , and perform the following recursive updates once at the end of of the overall algorithm for (P1) (i.e., once per minibatch of frames): Here C t denotes the final sparse codes estimates from the algorithm for (P1). In this case, when performing the inner update (9), the contributions of the two terms on the right hand side of (10) are incorporated separately. The matrices Q t ∈ C n×m and G t ∈ C m×m are small, constant-sized matrices 2 In our experiments, we set W to be the reshaped first column of the n×n identity matrix, which worked well. whose dimensions are independent of the time index t and the dimensions of the frame sequence, so they are stored for efficient use in the next minibatch. Moreover, the matrix C t is sparse, so all the matrix-matrix multiplications in (11) (or the matrix-vector multiplications in (10)) are computed using efficient sparse operations.

B. The Unitary Dictionary Variation
In the case of (P2), unlike for (P1), we do not perform block coordinate descent over the columns of D and (Z t ) H . Instead we minimize (P2) with respect to each of the matrices Z t and D directly and exploit simple closed-form solutions for the matrix-valued updates. The following subsections describe the solutions to the Z t and D update subproblems.

1) Sparse Coding
Step: Since D is a unitary matrix, minimizing (P2) with respect to Z t yields the following subproblem: where P t is again the matrix whose lth column is P l x t . The solution to (12) is given by the simple elementwise hardthresholding operation 2) Dictionary Update Step: Minimizing (P2) with respect to D results in the following optimization: where C j (Z j ) H is used for notational convenience. Using the definitions of the matricesP 1:t andC 1:t from the dictionary atom updates in Section IV-A, we can equivalently write (14) as Problem (15) is a well-known orthogonal Procrustes problem [48]. The solution is given bŷ where U ΣV H is a full SVD ofP 1:tC 1:t . Similarly as in (11), let Q t =P 1:tC 1:t . The matrix Q t ∈ C n×n can be recursively updated over time according to (11) (or (10)), so the dictionary update step can be performed efficiently and fully online.

C. Image Update Step
Minimizing (P1) or (P2) with respect to the minibatch x t yields the following quadratic sub-problem: Problem (17) is a least squares problem with normal equation In applications such as video denoising or inpainting, the matrix pre-multiplying x t in (18) is a diagonal matrix that can be efficiently pre-computed and inverted. More generally, in inverse problems where the matrix pre-multiplying x t (i.e., (18) is not diagonal or readily diagonalizable (e.g., in dynamic MRI with multiple coils), we minimize (17) by applying an iterative optimization method. One can use any classical algorithm in this case, such as the conjugate gradient (CG) method. For the experiments shown in Section V, we used a few iterations (indexed by i) of the simple proximal gradient method [49], [50] with updates of the form where θ(  (19) corresponds to a simple least squares problem that is solved efficiently by inverting a diagonal matrix I + 2τ i λ S M l=1 P T l P l (arising from the normal equation of the proximal operation), which is pre-computed. A constant step-size τ k = τ < 2/ A t 2 2 suffices for convergence [49]. Moreover, the iterations (19) monotonically decrease the objective in (17) when a constant step size τ ≤ 1/ A t 2 2 is used [50]. Fig. 2 summarizes the overall OnAIR algorithms for the online optimization Problems (P1) and (P2). In practice, we typically use more iterations K for reconstructing the first minibatch of frames, to create a good warm start for further efficient online optimization. The initialD 0 in the algorithm can be set to an analytical dictionary (e.g., based on the DCT or wavelets) and we set the initialẐ 0 to a zero matrix.

D. Computational Cost and Convergence
The computational cost for each time index t of the proposed algorithms for solving the online image reconstruction problems (P1) and (P2) scales as O(n 2 M ), where D ∈ C n×m with m ∝ n assumed, and M is the number of (overlapping) patches in each temporal window. The cost is dominated by various matrix-vector multiplications. Assuming each window's lengthM n, the memory (storage) requirement for the proposed algorithm scales as O(nM ), which is the space required to store the image patches of x t when performing the updates for (P1) or (P2). Typically the minibatch sizeM is small (to allow better tracking of temporal dynamics), so the number and maximum temporal width of 3D patches in each window are also small, ensuring modest memory usage for the proposed online methods.
While the overall computational cost and memory requirements for each time index t are similar for the algorithms for (P1) and (P2), the simple matrix-valued forms of the dictionary and sparse code updates for (P2) result in practice, in a severalfold decrease in actual runtimes due to optimizations inherent to matrix-valued computations in modern linear algebra libraries.
The proposed algorithms involve either exact block coordinate descent updates or for example, proximal gradient OnAIR Algorithms for (P1) and (P2) Inputs : measurement sequence {y t }, weights λ S and λ Z , rank r, upper bound L, forgetting factor ρ, number of dictionary learning iterationsK, number of proximal gradient iterationsK, and the number of outer iterations per minibatch K. Outputs : reconstructed dynamic image sequence {f t }. Initialization: Initial (first) estimates {f t 0 } of new incoming (measured) frames, initialD 0 andẐ 0 at t = 0, Q 0 = 0 n×m and G 0 = 0 n×m . For t = 1, 2, 3, . . ., do Warm Start: is set as follows: frames estimated in previous minibatches are set to a ρ-weighted average of those estimates and new frames are set as in {g t 0 }. For k = 1 : K repeat 1) Dictionary Learning Step: For (P1), update (D t k , C t k ) by sequentially updating the 1 ≤ i ≤ m coefficient and dictionary columns using (6) and (8) (that uses (9) and (10)) respectively, withK iterations and initial D t k−1 , C t k−1 . For (P2), update (D t k , C t k ) usingK alternations between performing (13) and (16) with initial D t k−1 , C t k−1 .

2) Image Update
Step: Update x t k using (18) or with K iterations of an optimization scheme (e.g., CG or (19)), with initialization x t k−1 . End Recursive Updates: Update Q t and G t using (11). Output Updates: SetẐ t = (C t K ) H ,D t = D t K andx t = x t K . For frames not occurring in future windows,f t is set to a ρ-weighted average of estimates over minibatches. End iterations (with appropriate step size) when the matrix premultiplying x t in (18) is not readily diagonalizable. These updates are guaranteed to monotonically decrease the objectives in (P1) and (P2) for each time index t. Whether the overall iterate sequence produced by the algorithms also converges over time (see e.g., [14]) is an interesting open question that we leave for future work.

V. NUMERICAL EXPERIMENTS
This section presents extensive numerical experiments illustrating the usefulness of the proposed OnAIR methods. We consider two inverse problem applications in this work: video reconstruction (inpainting) from noisy and subsampled pixels, and dynamic MRI reconstruction from highly undersampled k-t space data. The algorithm for (P1) is dubbed OnAIR-LD when the parameter r is lower than its maximum setting (i.e., low-rank dictionary atoms) and is dubbed OnAIR-FD (i.e., full-rank dictionary atoms) otherwise. The algorithm for (P2) is dubbed OnAIR-UD (unitary dictionary).   I: PSNR values in decibels (dB) for video inpainting for three videos from the BM4D dataset at various percentages of missing pixels. The methods considered are the proposed OnAIR-FD method (i.e., r = 5), the proposed OnAIR-UD method, online inpainting with a fixed DCT dictionary, the batch dictionary learning-based reconstruction method with r = 5, 2D (frame-by-frame cubic) interpolation, and 3D interpolation. The best PSNR value for each undersampling level for each video is in bold.
subsections discuss the applications and results. A link to software to reproduce our results will be provided at http: //web.eecs.umich.edu/ ∼ fessler/.
A. Video Inpainting 1) Framework: First, we consider video inpainting or reconstruction from subsampled and potentially noisy pixels. We work with the publicly available videos 3 provided by the authors of the BM4D method [51]. We process the first 150 frames of each video at native resolution. We measure a (uniform) random subset of the pixels in each frame of the videos, and also simulate additive (zero mean) Gaussian noise for the measured pixels in some experiments. The proposed OnAIR methods are then used to reconstruct the videos from their corrupted (noisy and/or subsampled) measurements.
The parameters for the proposed OnAIR methods are chosen as follows. We used a sliding (temporal) window of length M = 5 frames with a temporal stride of 1 frame, to reconstruct (maximally overlapping) minibatches of frames. In each window, we extracted 8 × 8 × 5 overlapping spatiotemporal patches with a spatial stride of 2 pixels. We learned a square 320×320 dictionary, and the operator R(·) reshaped dictionary atoms into 64 × 5 space-time matrices. We ran the OnAIR algorithms for (P1) and (P2) for K = 7 iterations in each temporal window (minibatch), withK = 1 inner iteration (of block coordinate descent or alternation) for updating (D, Z t ), and used a direct reconstruction as per (18) in the image update step. A forgetting factor ρ = 0.9 was observed to work well. We ran the algorithms for more (K = 50) iterations for the first minibatch of data, for which the initial D was the discrete cosine transform (DCT) matrix and the initial sparse codes were zero. Newly arrived frames were first initialized (i.e., in the first minibatch they occur) using 2D cubic interpolation. We simulated various levels of subsampling of the videos (with and without noise), and we chose r = 5 or full-rank atoms in (P1), which outperformed low-rank atoms in the experiments here. The videos in this section have substantial temporal motion, so allowing full-rank atoms enabled the algorithm to better learn the dynamic temporal features of the data. Section V-B demonstrates the usefulness of low-rank atoms in (P1). We tuned the weights λ S and λ Z for (P1) and (P2) to achieve good reconstruction quality at an intermediate 3 The data is available at http://www.cs.tut.fi/ ∼ foi/GCF-BM3D/.
undersampling factor (70% missing pixels) for each video, and used the same parameters for other factors.
We measure the performance of our methods using the (3D) peak signal-to-noise ratio (PSNR) metric that is computed as the ratio of the peak pixel intensity in the video to the root mean square reconstruction error relative to the ground truth video. All PSNR values are expressed in decibels (dB).
We compare the performance of the proposed OnAIR-FD (r = 5) and OnAIR-UD methods for (P1) and (P2) with that of an identical online algorithm but which uses a fixed DCT dictionary. We also produce reconstructions using a "batch" version of the method for (P1) [13] (with r = 5) that processes all video frames jointly. 4 For each comparison method, we used the same patch dimensions, initializations, etc., and tuned the parameters λ S and λ Z of each method individually. The one exception is that we used a spatial patch stride of 4 pixels for the Batch Learning method rather than the 2 pixel spatial stride used for the online methods. This change was made because the batch method is memory intensiveit extracts and processes image patches from the entire video concurrently-so it was necessary to process fewer patches to make the computation feasible. We ran the batch learningbased reconstruction for 20 iterations, which worked well. Finally, we also compare the OnAIR reconstructions to baseline reconstructions produced by 2D (frame-by-frame cubic) interpolation and 3D interpolation (using the natural neighbor method in MATLAB).
2) Results: Table I lists the PSNR values (in dB) for the various reconstruction methods for different levels of subsampling (from 50% to 90% missing pixels) of three videos (without added noise) from the BM4D dataset. Table II shows analogous results for the Coastguard video when i.i.d. zero-mean Gaussian noise with 25 dB PSNR was added before sampling the pixels. The proposed OnAIR-FD method typically provides the best PSNRs at higher undersampling factors, and the proposed OnAIR-UD method (with a more structured/unitary 5 D) performs better at lower undersampling factors.  The best PSNRs achieved by the OnAIR methods are typically better (by up to 1.4 dB) than for the Batch Learning scheme. Both OnAIR variations also typically outperform (by up to 2 dB) the online method with a fixed DCT dictionary (the initial D in our methods), suggesting that the dictionaries learned by the proposed methods are successfully adapting to and reconstructing underlying features of the videos.
Figs. 3 and 4 show the original and reconstructed frames for a few representative frames of the Coastguard and Flower Garden videos. Results for multiple methods are shown. Fig. 3 shows that the proposed OnAIR-FD method produces visually more accurate reconstructions of the texture in the waves and   5) and OnAIR-UD schemes, online inpainting with a fixed DCT dictionary, batch dictionary learning-based reconstruction (r = 5), 2D (frame-by-frame cubic) interpolation, and 3D interpolation. The best PSNR for each undersampling is in bold. , OnAIR-UD, the online method with fixed DCT dictionary, the batch learningbased method, and 3D interpolation, for the Coastguard video with 70% missing pixels. also produces fewer artifacts near the boats in the water and the rocks on the shore. Fig. 4 illustrates that the proposed OnAIR method produces a sharper reconstruction with less smoothing artifacts than the online method with a fixed DCT dictionary and the batch learning-based reconstruction method, and it is less noisy than the interpolation-based reconstruction.
3) Properties: Fig. 5 shows the frame-by-frame (2D) PSNRs for the Coastguard video inpainted from 70% missing pixels, using various methods. Clearly the proposed OnAIR-FD method achieves generally higher PSNRs across frames of the video. The overall trends in PSNR over frames are similar across the methods and are due to motion in the original videos, with more motion generally resulting in lower PSNRs. Fig. 6 shows a representative example of a (final) learned dictionary produced by the proposed OnAIR-FD method for the Bus video, along with the initial (at t = 0) DCT dictionary. The dictionaries each contain 320 atoms, each of which is an 8 × 8 × 5 space-time tensor. We visualize each atom by plotting the first 8 × 8 x − y slice of the atom and also plotting the y − t profile from a vertical slice through the middle of each atom tensor. The x − y (first slice) images show that the learned dictionary has adapted to both smooth and sharp gradients in the image, and the dynamic (evolved) nature of the y −t profiles shows that the dictionary atoms have adapted Initial (1st slice) Initial (y-t) Learned (1st slice) Learned (y-t) to temporal trends in the data.
To assess the relative efficiency of the OnAIR methods compared to the Batch Learning method, we measured runtimes on the Coastguard video with the same fixed patch sizes, patch strides, and numbers of iterations for each method. 6 In this case, the OnAIR-FD and OnAIR-UD methods were 1.8x and 4.0x faster, respectively, than the Batch Learning method. In addition, the OnAIR methods typically require much fewer iterations (per minibatch) compared to the batch method to achieve their best reconstruction accuracies, so, in practice, when utilizing half as many outer iterations, the OnAIR-FD and OnAIR-UD methods were approximately 3.6x and 8.0x faster, respectively, than the batch scheme.
B. Dynamic MRI Reconstruction 1) Framework: Here, we demonstrate the usefulness of the proposed OnAIR methods for reconstructing dynamic MRI data from highly undersampled k-t space measurements. We work with the multi-coil (12-element coil array) cardiac perfusion data [24] and the PINCAT data [12], [52] from prior works. For the cardiac perfusion data, we retrospectively undersampled the k-t space using variable-density random Cartesian undersampling with a different undersampling pattern for each time frame, and for the PINCAT data we used pseudo-radial sampling with a random rotation of radial lines between frames. The undersampling factors tested are shown in Table III. For each dataset, we obtained reconstructions using the proposed OnAIR methods, the online method but with a fixed DCT dictionary, and the batch learning-based reconstruction (based on (P1)) method [13]. We also ran the recent L+S [24] and k-t SLR [53] methods, two sophisticated batch methods for dynamic MRI reconstruction that process all the frames jointly. Finally, we also computed a baseline reconstruction in each experiment by performing zeroth order interpolation across time at non-sampled k-t space locations (by filling such locations with the nearest non-zero entry along time) and then backpropagating the filled k-t space to image space by pre-multiplying with the A H corresponding to fully sampled data. The first estimates (g t 0 ) of newly arrived frames in our OnAIR methods are also computed via such zeroth order interpolation, but using only the estimated (i.e., already once processed and reconstructed) nearest (older) frames.
For the online schemes, we used 8 × 8 × 5 spatiotemporal patches withM = 5 frames per temporal window (minibatch) and a (temporal) window stride of 1 frame. We extracted overlapping patches in each minibatch using a spatial stride of 2 pixels along each dimension. We learned a square 320 × 320 dictionary whose atoms when reshaped into 64 × 5 spacetime matrices had rank r = 1, which worked well in our experiments. A forgetting factor of ρ = 0.9 was observed to work well. We ran the online schemes for K = 7 outer iterations 7 per minibatch, withK = 1 iteration in the dictionary learning step andK = 10 proximal gradient steps in the image update step, respectively. We used K = 50 outer iterations for the first minibatch to warm start the algorithm and used an initial DCTD 0 , and the initial sparse codes were zero. After one complete pass of the online schemes over all the frames, we performed another pass over the frames, using the reconstructed frames and learned dictionary from the first pass as the first initializations in the second pass. The additional pass gives a slight image quality improvement over the first (a single pass corresponds to a fully online scheme) as will be shown later.
For the batch dictionary learning-based reconstruction method, we used the same patch dimensions, strides, and initializations as for the online methods. We ran the batch method for 50 iterations. For the L+S and k-t SLR methods, we used the publicly available MATLAB implementations from [54] and [55], respectively, and ran each method to convergence. The regularization parameters (weights) for all the methods here were tuned for each dataset by sweeping them over a range of values and selecting values that achieved good reconstruction quality at intermediate k-t space undersampling factors. We measured the dMRI reconstruction quality using the normalized root mean square error (NRMSE) metric expressed as percentage that is computed as wherex is a candidate reconstruction and x ref is the reference reconstruction (e.g., computed from "fully" sampled data).
2) Results and Comparisons: Table III shows the reconstruction NRMSE values obtained using various methods for the cardiac perfusion and PINCAT datasets at several undersampling factors. The proposed OnAIR-LD method achieves lower NRMSE values in almost every case compared to the L+S and k-t SLR methods, the online DCT scheme, the baseline reconstruction, and the batch learning-based reconstruction scheme. Unlike the batch schemes (i.e., L+S, k-t SLR, and the batch learning-based reconstruction) that process or learn over all the data jointly, the OnAIR-LD scheme only directly processes data corresponding to a minibatch of 5 frames at any time. Yet OnAIR-LD outperforms the other methods because of its memory (via the forgetting factor) and local temporal adaptivity (or tracking). These results show that the proposed OnAIR methods are well-suited for processing streaming data. Fig. 7 shows the reconstructions and reconstruction error maps (magnitudes displayed) for some representative frames from the cardiac perfusion and PINCAT datasets at 12x and 7x undersampling, respectively. The error maps indicate that the proposed OnAIR-LD scheme often produces fewer artifacts compared to the existing methods. Fig. 8 shows that the OnAIR-LD scheme also typically achieves better frame-by-frame NRMSE compared to the other methods. Finally, Fig. 9 shows the reconstructed y − t profiles for various methods obtained by extracting the same vertical line segment from each reconstructed frame of the PINCAT data and concatenating them. The online DCT scheme and the batch methods L+S and k-t SLR show line-like or additional smoothing artifacts that are not produced by the proposed OnAIR method.
3) Properties: Table IV investigates the properties of the proposed OnAIR methods in more detail using the cardiac perfusion data. Specifically, it compares the NRMSE values produced by the OnAIR-LD scheme with one or two passes over the data, the OnAIR-UD method, and the online method with a fixed DCT dictionary. In addition, we ran the online method but with a fixed "oracle" dictionary learned from patches of the reference (true) reconstruction by solving the DINO-KAT learning problem (2). The oracle dictionary was computed based on the "fully" sampled data, so it can be viewed as the "best" dictionary that one could learn from the undersampled dataset. From Table IV, we see that the NRMSE values achieved by the OnAIR-LD scheme with two passes are within 0.0% -0.5% of the oracle NRMSE values, which suggests that the proposed scheme is able to learn dictionaries with good representational (recovery) qualities from highly undersampled data. Moreover, the performance of the OnAIR-LD scheme with a single pass is almost identical to that with two passes, suggesting promise for fully online dMRI reconstruction. The OnAIR-LD method outperformed the OnAIR-UD scheme for the cardiac perfusion data indicating that temporal low-rank properties better characterize the dataset than unitary models. All the OnAIR schemes outperformed the nonadaptive yet online DCT-based scheme.
To assess the relative efficiency of the OnAIR methods com-  TABLE III: NRMSE values as percentages at several undersampling factors for the cardiac perfusion data with Cartesian sampling (left) and for the PINCAT data with pseudo-radial sampling (right). The methods compared are the proposed OnAIR-LD (r = 1) method, the online scheme with a fixed DCT dictionary, the batch learning-based reconstruction (r = 1), the L+S method, the k-t SLR method, and a baseline reconstruction. The best NRMSE achieved for each undersampling and each dataset is in bold. The right four columns of each panel depict the corresponding reconstruction error magnitudes (w.r.t. reference) for OnAIR-LD, the online method with a fixed DCT dictionary, the k-t SLR method, and the L+S method, respectively. Compared to the competing methods, OnAIR-LD achieves (3D) NRMSE improvements of 0.6 dB, 1.9 dB, and 0.7 dB, respectively, and 0.4 dB, 1.0 dB, and 2.3 dB, respectively, for the cardiac perfusion and PINCAT data. , the online method with fixed DCT dictionary, and the k-t SLR and L+S methods, for the PINCAT data with 9x undersampling and pseudo-radial sampling.
pared to the Batch Learning method, we measured runtimes on the PINCAT dataset with the same parameter settings for each method used to generate the results in Table III. 8 With these settings, the OnAIR-LD and OnAIR-UD methods were 9.8x and 40.4x faster, respectively, than the Batch Learning method. One key reason for these significant performance improvements is that 50 outer iterations were required by 8 Experiments were performed in MATLAB R2017b on a 2016 MacBook Pro with a 2.7 GHz Intel Core i7 processor and 16 GB RAM.

Reference
OnAIR-LD Online (DCT) k-t SLR L+S Fig. 9: Temporal (y − t) profiles of a spatial vertical line cross section for the reference PINCAT reconstruction, the proposed OnAIR-LD (r = 1) method, the online method with fixed DCT dictionary, the k-t SLR method, and the L+S method for 14x undersampling (pseudo-radial sampling).
the batch method while only 7 outer iterations were required by the online methods (after the first minibatch) to achieve the reported accuracies; this result suggests that the online dictionary adaptation performed by the OnAIR methods is both computationally efficient and allows the model to better adapt to the underlying structure of the data. The additional speedup of the OnAIR-UD method with respect to OnAIR-LD is attributable to the relative efficiency of the matrix-valued (D, B) updates of the unitary method compared to the less optimized block coordinate descent iterations over the columns  The methods compared are an oracle online scheme (where the dictionary is learned from patches of the reference and fixed during reconstruction), the proposed OnAIR-LD (r = 1) method with a single pass or two passes over the frames, the proposed OnAIR-UD method, and the online scheme with a fixed DCT dictionary. The best NRMSE for each undersampling is in bold.
Initial Learned (real) Learned (imaginary) Fig. 10: Dictionaries for the PINCAT data with 9x undersampling. Left: the atoms of the initial DCT dictionary. Center and Right: the real and imaginary parts of the final learned dictionary produced by OnAIR-LD (r = 1). The dictionary atoms are 8 × 8 × 5 tensors, and only the first 8 × 8 (x − y) slice of each atom is displayed.
of D and B prescribed by OnAIR-LD. Fig. 10 shows an example of a learned dictionary produced by the proposed OnAIR-LD method on the PINCAT dataset at 9x undersampling, which is compared with the initial DCT dictionary. The dictionaries have 320 atoms, each a 8 × 8 × 5 complex-valued space-time tensor. Since the OnAIR LD dictionary atoms have rank r = 1 when reshaped into 64 × 5 space-time matrices, we directly display the real and imaginary parts of the first 8 × 8 (x − y) slice of each learned atom. The initial DCT is also displayed similarly. The (eventual) learned dictionary has clearly evolved significantly from the initial DCT atoms and has adapted to certain smooth and sharp textures at various orientations in the data.
Recall that we chose full-rank (r = 5) atoms in the video inpainting experiments, while here we chose low-rank (r = 1) atoms. Intuitively, low-rank atoms are a better model for the dynamic MRI data because the videos have high temporal correlation and rank-1 atoms are necessarily constant in their temporal dimension. Conversely, the videos from Section V-A contained significant camera motion and thus dictionary atoms with more temporal variation (i.e., higher rank) enabled more accurate reconstructions.

VI. CONCLUSIONS
This paper has presented a framework for online estimation of dynamic image sequences by learning dictionaries. Various properties were also studied for the learned dictionary such as a unitary property, and low-rank atoms, which offer additional efficiency or robustness to artifacts. The proposed OnAIR algorithms sequentially and efficiently update the images, dictionary, and sparse coefficients of image patches from streaming measurements. Importantly, our algorithms can process arbitrarily long video sequences with constant memory usage over time. Our numerical experiments demonstrated that the proposed methods produce accurate reconstructions for video inpainting and dynamic MRI reconstruction. The proposed methods may also be suitable for other inverse problems, including medical imaging applications such as interventional imaging, and other video-processing tasks from computer vision. We hope to investigate these application domains as well as demonstrate potential real-time applicability for OnAIR approaches in future work.