Compressive Imaging by Generalized Total Variation Minimization Jie - PowerPoint PPT Presentation

Compressive Imaging by Generalized Total Variation Minimization Jie Yan and Wu-Sheng Lu Department of Electrical and Computer Engineering University of Victoria, Victoria, BC, Canada Aug 10, 2012 1 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 2 / 26

Introduction Signal recovery from compressive sensing (CS) measurements conventionally achieved by minimizing the ℓ 1 norm. Several authors recently reported nonconvex ℓ p regularization improves recovery performance, using a p less than 1. Iteratively reweighted ℓ 1 minimization (IRL1) proves successful in tackling the ℓ p problem without having to perform nonconvex optimization. 3 / 26

Introduction A conventional total variation (TV) minimization for image recovery TV ( U ) s.t. � R ◦ ( F U ) − B � 2 F < σ 2 min (1) U where ◦ is the Hadamard product operator, F denotes the 2-D Fourier transform operator, R represents a random sampling matrix whose entries are either 1 or 0, and B stores the compressive sampled measurements. Inspired by the relationship between ℓ p and ℓ 1 norms, we generalize the concept of TV to a p th-power type TV with 0 ≤ p ≤ 1 . Inspired by the IRL1, we propose an iteratively reweighted TV minimization algorithm to approach the solution that minimizes the p th-power TV. 4 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 5 / 26

Anisotropic TV Anisotropic total variation (TV) of a digital image U ∈ R n × n is defined as n − 1 n − 1 � � TV ( U ) = ( | U i , j − U i + 1 , j | + | U i , j − U i , j + 1 | ) i = 1 j = 1 (2) n − 1 n − 1 � � + | U i , n − U i + 1 , n | + | U n , j − U n , j + 1 | i = 1 j = 1 Under the periodic boundary condition, TV can be expressed in the form TV ( U ) = � DU � 1 + � UD T � 1 (3) with D ∈ R n × n as a circulant matrix with the first row [ 1 − 1 0 · · · 0 ] . The notation � X � 1 denotes the sum of magnitudes of all the entries in X , i.e., � | x i , j | . 6 / 26

Generalized p th power TV We extend the concept of TV by defining a generalized p th power TV (GTV) as TV p ( U ) = � DU � p + � UD T � p The newly introduced notation � X � p resembles an ℓ p norm as it expresses the sum of p th power magnitudes of all the entries in X , i.e., � | x i , j | p . We are inspired to investigate the GTV, i.e., TV p ( U ) , with 0 ≤ p ≤ 1 , for CS image recovery. 7 / 26

Weighted TV However attractive TV p appears in promoting a sparser TV, optimization of TV p -related problem is nonconvex, nonsmooth, and so far no algorithms have been discovered in handling such a problem to the authors’ knowledge. Instead, we introduce a weighted TV (WTV) as TV w ( U ) = � W x ◦ ( DU ) � 1 + � W y ◦ ( UD T ) � 1 TV w ( U ) becomes the conventional TV ( U ) when all entries in W x and W y are ones. 8 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 9 / 26

A Power-Iterative Reweighting Strategy We adopt a power-iterative strategy to approach the solution for a TV 0 -regularized problem. Set p = 1 , l = 1 , W x = W y = 1 . 1 Solve the WTV-regularized problem for U ( l ) 2 TV w ( U ) s.t. � R ◦ ( F U ) − B � 2 F < σ 2 min (4) U Terminate if p = 0 ; otherwise, set p = p − 0 . 1 and update the 3 weights W x and W y as W x = | DU ( l ) + ǫ | . p − 1 , W y = | U ( l ) D T + ǫ | . p − 1 (5) Then set l = l + 1 and repeat from step 2. 10 / 26

A Power-Iterative Reweighting Strategy Note that by Eq. (5), TV w ( U ) essentially becomes TV p ( U ) for U in a neighborhood of iterate U ( l ) . Consequently, nonconvex minimization of TV p ( U ) can practically be achieved by a series of convex minimization of TV w ( U ) using the above strategy. 11 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 12 / 26

WTV-Regularized Minimization The analysis has led to the WTV-regularized problem TV w ( U ) s.t. � R ◦ ( F U ) − B � 2 F < σ 2 min (6) U We propose to solve the problem using a Split Bregman approach, but with important changes. Using Bregman iteration, we reduce the problem to TV w ( U ) + µ U ( k + 1 ) = argmin 2 � R ◦ ( F U ) − B ( k ) � 2 (7a) F U B ( k + 1 ) = B ( k ) + B − R ◦ ( F U ( k + 1 ) ) (7b) 13 / 26

WTV-Regularized Minimization A splitting strategy applied to (7a) leads to the formulation U � W x ◦ D x � 1 + � W y ◦ D y � 1 + µ 2 � R ◦ ( F U ) − B ( k ) � 2 min (8a) F s.t. D x = DV , D y = VD T , U = V (8b) Here we apply the Split Bregman to split D x = DV and D y = VD T , but also introduce an additional split as U = V . Such a split allows us to decompose the most computationally expensive step of the algorithm into two much simpler steps. 14 / 26

WTV-Regularized Minimization Enforcing the constraints in (8b), we mimize the following function with respect to { U , V , D x , D y } min � W x ◦ D x � 1 + � W y ◦ D y � 1 + µ 2 � R ◦ ( F U ) − B ( k ) � 2 F + λ F + λ 2 � D y − VD T − E ( h ) 2 � D x − DV − E ( h ) x � 2 y � 2 F + ν 2 � U − V − G ( h ) � 2 F E ( h ) x , E ( h ) and G ( h ) are updated through Bregman iterations. y 15 / 26

WTV-Regularized Minimization At the h th iteration we solve four subproblems µ U ( h + 1 ) = argmin 2 � R ◦ ( F U ) − B ( k ) � 2 F U (9) + ν 2 � U − V ( h ) − G ( h ) � 2 F ν V ( h + 1 ) = argmin 2 � V − U ( h + 1 ) + G ( h ) � 2 F V (10) + λ F + λ 2 � VD T + E ( h ) 2 � DV + E ( h ) − D ( h ) − D ( h ) x � 2 y � 2 x y F � W x ◦ D x � 1 + λ 2 � D x − DV ( h + 1 ) − E ( h ) D ( h + 1 ) x � 2 = argmin (11a) x F D x � W y ◦ D y � 1 + λ 2 � D y − V ( h + 1 ) D T − E ( h ) D ( h + 1 ) y � 2 = argmin (11b) y F D y 16 / 26

WTV-Regularized Minimization The problems in (11a) and (11b) can be solved simply by soft shrinkage operations as = T W x /λ ( DV ( h + 1 ) + E ( h ) D ( h + 1 ) x ) (12a) x = T W y /λ ( V ( h + 1 ) D T + E ( h ) D ( h + 1 ) y ) (12b) y Solving problems (9) and (10) are however far from trivial 17 / 26

WTV-Regularized Minimization We write first-order optimality condition of (9) as µ F T R ◦ F U + ν U = µ F T R ◦ B k + ν ( V ( h ) + G ( h ) ) (13) Multiply both sides of (13) by F on the left, we obtain solution of (13) as U ( h + 1 ) = F T �� µ R ◦ B k + ν F ( V ( h ) + G ( h ) ) � � ◦ / ( µ R + ν ) (14) 18 / 26

WTV-Regularized Minimization To solve (10), we write its first-order optimality condition as ν V + λ D T DV + λ VD T D = C ( h ) (15) where C ( h ) = ν ( U ( h + 1 ) − G ( h ) ) (16) + λ D T ( D ( h ) − E ( h ) x ) + λ ( D ( h ) − E ( h ) y ) D x y Circulant matrix D can be diagonalized by the 2-D Fourier transform F as D = F T Λ F . 19 / 26

WTV-Regularized Minimization Multiply both sides of (15) by F on the left and F T on the right, we have ν ˜ V + λ ( T ˜ V + ˜ VT ) = F C ( h ) F T (17) V = F V F T and T = Λ ∗ Λ . where ˜ T is a diagonal matrix, we can further express (17) as ( ν + λ T r + λ T c ) ◦ ˜ V = F C ( h ) F T (18) where T r has each element in its i th row as T i , i while T c has each element in its i th column as T i , i . Therefore, we obtain solution of (10) as V ( h + 1 ) = F T � � ( F C ( h ) F T ) ◦ / ( ν + λ T r + λ T c ) F (19) 20 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 21 / 26

MRI of the Shepp-Logan Phantom A normalized Shepp-Logan phantom of size 256 × 256 , was measured at 2521 locations in the 2D Fourier plane ( k -space). The sampling pattern was a star-shaped pattern consisting of only 10 radial lines. Recover the image based on the 2521 star-shaped 2D Fourier samples. 22 / 26

MRI of the Shepp-Logan Phantom (a) Star-shaped sampling pattern (b) Minimum energy reconstruction (c) Minimum TV reconstruction (d) Minimum GTV reconstruction with p = 0 23 / 26

Recovered Shepp-Logan Phantom Minimum GTV reconstruction with p = 0 The signal to noise (SNR) ratio: 16.3 dB. Computation time on a PC laptop with a 2.67 GHz Intel quad-core processor: 770.7 seconds. Minimum TV reconstruction SNR: 8.8 dB. Computation time: 756.8 seconds. 24 / 26

OUTLINE Introduction 1 Generalized TV and Weighted TV 2 A Power-Iterative Reweighting Strategy 3 WTV-Regularized Minimization: The Splitting Technique 4 5 Performance on Image Reconstruction Conclusions 6 25 / 26