Distributed Reconstruction via Alternating Direction Method

With the development of compressive sensing theory, image reconstruction from few-view projections has received considerable research attentions in the field of computed tomography (CT). Total-variation- (TV-) based CT image reconstruction has been shown to be experimentally capable of producing accurate reconstructions from sparse-view data. In this study, a distributed reconstruction algorithm based on TV minimization has been developed. This algorithm is very simple as it uses the alternating direction method. The proposed method can accelerate the alternating direction total variation minimization (ADTVM) algorithm without losing accuracy.


Introduction
Image reconstruction algorithms implemented in existing computed tomography (CT) scanners require projection data to be available in proportional space [1,2]. However, in CT imaging of biological specimens, data collection at a large number of projection views may result in radiationinduced object deformation. Recently, methods based on the corresponding constrained total variation (TV) or 1 -norm minimization have been widely studied in reconstruction from sparse-view data [3][4][5][6][7][8][9]. Sidky et al. [5] proposed the adaptive steepest descent projection onto convex sets (ASD-POCS) algorithm for CT image reconstruction. This method can realize exact-image reconstruction using fewer measurements. The alternating direction method (ADM) is an efficient approach for optimization problems. And the split Bregman-TV method based on alternating Bregman iterative approach was proposed and converged well in [8] as a solution for sparse-view CT reconstruction. A TV minimization iterative algorithm using the ADM based on augmented Lagrangian function was also proposed [10,11]. Li et al. proposed a more robust and efficient algorithm nonmonotone alternating direction algorithm (NADA) [12] in 2012, which integrates alternating direction and nonmonotone line search. An alternating direction total variation minimization (ADTVM) algorithm for few-views reconstruction [13] was developed inspired by the literature [10][11][12]. The augmented Lagrangian function-based ADM is actually equivalent to the Bregman iterative method when the constraints are linear [14]. However, the expression in [13] is simpler than that in [8].
The CT image reconstruction problem is a large-scale problem. The ADTVM algorithm [13] is not directly suitable for distributed implementation. Boyd et al. [15] argued that the alternating direction method of multipliers is well suited for distributed convex optimization, in particular, for largescale problems arising in statistics, machine learning, and other related areas. In this study, a distributed algorithm called distributed alternating direction total variation minimization (Dis-ADTVM) is developed using ADM. The proposed algorithm is as simple as the ADTVM algorithm and can accelerate the latter without accuracy loss.

Method
The imaging model can be approximated using the following discrete linear system: where the vector has length , which is the number of measured projection rays; the vector ⃗ has length , which is the number of expansion elements used in representing the object function ( ⃗ ); the system matrix is a pixel-driven projection operator. Sparse-views projection data are not sufficient for exact reconstruction. The problem we consider in this study is illposed. To solve linear system (1), we use a regularization method with anisotropic TV minimization, as follows: where ‖ ⃗ ‖ ≜ ∑ ‖ ⃗ ‖ 1 and denotes the differential operator along direction . In particular, 1 and 2 denote the horizontal and vertical differential operators, respectively, for two-dimensional form. and are separated in along vertical direction as min ⃗ We consider a variant of (3) as follows: where ⃗ denotes ⃗ in node and ⃗ = ⃗ , = 1, 2, . . . , , in node . Its corresponding augmented Lagrangian function is where is Lagrange multiplier and the parameters and are both used to balance the terms. The ADM is used to  solve the problem that minimizes the augmented Lagrangian function; that is, The final algorithm of Dis-ADTVM can be expressed as where ⃗ +1 and +1 can be computed in node and + stands for the Moore-Penrose pseudoinverse of matrix .
Computing the pseudoinverse at each iteration is too costly to implement numerically, while we use NADA [12] to solve " -subproblem" in (7). The convergence analysis of ADTVM algorithm has been well analyzed in the literature [12], and the convergence proof of the distributed algorithms based on ADM can be found in [15,16]. The iterative algorithms using the ADM based on augmented Lagrangian function decompose the optimization problem into some simple subproblems with closed form solution. Therefore, the algorithms are efficient and practical for the low cost in each iteration. The NADA algorithm enables taking full advantages of the low-cost minimization in "easy" direction and allows relatively quick and large steps in the "hard" direction. The distributed algorithms can distribute some computation to individual nodes; thereby, the algorithms reduce the running time through data distribution and computation. The proposed Dis-ADTVM algorithm in this paper integrates above advantages and its derivation and implementation are as simple as the ADTVM algorithm.

Simulation Studies.
We perform numerical experiments to demonstrate and validate the proposed method for sparseview image reconstruction. Scanning and reconstruction parameters are listed in Table 1. Detector elements are equidistantly spaced at 0.127 mm.   The proposed method is compared with ASD-POCS algorithm [6] and the ADTVM algorithm [13], using the same parameters to validate their performance. Dual core is used in implementing the proposed distributed algorithm on two nodes.
In the experimental configuration, we use one detector for data acquisition by taking 36 angular samples evenly distributed over an angular range of 360 ∘ . The size of the phantom simulation is set as follows. Image size is 256 × 256 = 65536 voxels, and projection data size is 36 views, with 512 detectors or 18432 measured rays.
The images reconstructed from this set of data using the ASD-POCS algorithm, ADTVM algorithm, and the proposed distributed algorithm are shown in Figure 1. The profiles of these images along the central horizontal and vertical rows are presented in Figure 2. The number of iterations for the three algorithms is 1000 each. The parameters of ASD-POCS are same as those in [6]. The parameters  in the ADTVM algorithm and the proposed distributed algorithm are both = 1/4000 and = 32/4000.
We use the root mean squared error (RMSE) as a measure of the reconstruction error to demonstrate reconstruction accuracy quantitatively. The RMSE is defined as where and are the ideal phantom and the reconstruction, respectively, and is the total number of pixels in the image. The RMSEs of the reconstructions of the Shepp-Logan phantom are calculated. The results of the three methods are illustrated in Figure 3. Table 2 shows the RMSE of the reconstructions from the projection data above with the three algorithms. It is clear that the accuracy and effectivity of the ADTVM algorithm and the proposed distributed algorithm are both better than those of the ASD-POCS algorithm. This is due to the use of ADM and NADA algorithms. Moreover, we can see that the accuracies of the ADTVM algorithm and proposed distributed algorithm are both almost the same. This is because the derivation and implementation of distributed algorithm are very similar as the ADTVM algorithm.
The running time of the three algorithms is shown in Table 3 and Figure 4 for the phantom results in the aforementioned configuration. Timing is implemented based on the average of 10 computations. The speedup is approximately 1.4 on the average.

Reconstruction Using Real Data.
We perform experiments to reconstruct a head model from real data to further validate the proposed algorithm. Scanning and reconstruction parameters are listed in Table 4. Detector elements are equidistantly spaced at 0.635 mm.
We reconstruct a -axial slice for convenience. Images reconstructed using the ASD-POCS algorithm, the ADTVM algorithm, and the proposed distributed algorithm are shown in Figure 5. The numbers of iterations for the three algorithms are 200 and 1000, respectively. The experimental result suggests that ADTVM algorithm and the proposed distributed algorithm produce better reconstruction than ASD-POCS. Hence, the results of ADTVM algorithm and the proposed distributed algorithm are almost the same.
The running time of the three algorithms is shown in Table 5 and Figure 6 for the reconstruction results in the aforementioned configuration. Timing is implemented based on the average of 10 computations. The speedup is approximately 1.4 on the average. We use Amdahl's law [17] to predict the theoretical maximum speedup as follows: where denotes the fraction of the algorithm which is strictly serial. In the ADTVM algorithm, is about 0.1, so the theoretical maximum speedup in two nodes is = 1.82. The speedup in real experiments will be less than as the cost on communication exists in every iteration. The average speedup of the proposed distributed algorithm in all the experiments shows that the algorithm reduces the running time obviously.
All experiments are performed using C programming language under Visual Studio 2012 and OpenMP running on an AMAX Tesla workstation with Intel Xeon E5520 dual-core CPU 2.27 GHz and 24 GB memories. We partly refer to the MATLAB solver of "TVAL3" [11] for the implementation.

Conclusions
The Dis-ADTVM algorithm is as simple as the ADTVM algorithm and can accelerate the latter without accuracy loss. The new algorithm is well suited for CT image sparse-view reconstruction problem as a large-scale problem. It is clear that the Dis-ADTVM algorithm can be applied to other tomographic imaging modalities with linear system models. We will study the relationship between the performance and the number of nodes in a forthcoming paper.