DLCT-Toolbox, a Matlab package for the
Dictionary Learning Approach to Tomographic
Image Reconstruction
Sara Soltani ∗
Department of Applied Mathematics and Computer Science,
Technical University of Denmark, DK-2800 Kgs. Lyngby,
Denmark.
August 25, 2015
Contents
1 Introduction
1.1 Overall . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Toolbox and Software Dependencies
. . . . . . . . . . . . . . . .
1.3 Test Images . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.4 Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2
2
3
3
3
4
4
4
4
6
6
7
7
7
8
8
9
11
11
11
2 Package Details
2.1.1
2.1.2
2.1.3
2.1.4
2.1 Dictionary Learning Problem . . . . . . . . . . . . . . . . . . . .
DL_demo.m . . . . . . . . . . . . . . . . . . . . . . . . . .
DicLear_M.m . . . . . . . . . . . . . . . . . . . . . . . . .
Lea_Algorithm_D2.m . . . . . . . . . . . . . . . . . . . .
Lea_Algorithm_Dinf.m . . . . . . . . . . . . . . . . . . .
2.2 Mean Approximation Error . . . . . . . . . . . . . . . . . . . . .
MAEM_demo.m . . . . . . . . . . . . . . . . . . . . . . . . .
MAE_M.m . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 The Reconstruction Problem . . . . . . . . . . . . . . . . . . . .
RecM_demo.m . . . . . . . . . . . . . . . . . . . . . . . . .
RecM_Algorithm.m . . . . . . . . . . . . . . . . . . . . . .
2.4 The Tensor Dictionary Learning problem . . . . . . . . . . . . .
DL_Tensor_demo.m . . . . . . . . . . . . . . . . . . . . . .
DicLear_T.m . . . . . . . . . . . . . . . . . . . . . . . . .
2.3.1
2.3.2
2.4.1
2.4.2
2.2.1
2.2.2
∗ssol@dtu.dk, sarahsoltani@gmail.com
This work is part of the project HD-Tomo funded by Advanced Grant No. 291405 from the
European Research Council.
2010 Mathematics Subject Classification: 65F22, 65K10, 15A69.
Key words and Phrases: Tomography, Dictionary learning, Tensor decomposition, Inverse
problem, Regularization, Sparse representation, Image reconstruction.
1
S. Soltani
A Matlab Package for DL Approach to CT
2.4.3
2.4.4
Lea_Algorithm_T_D2.m . . . . . . . . . . . . . . . . . . .
Lea_Algorithm_T_Dinf.m . . . . . . . . . . . . . . . . . .
2.5 The Approximation Error by the Tensor Dictionary . . . . . . .
MAET_demo.m . . . . . . . . . . . . . . . . . . . . . . . . .
MAE_T.m . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.6 Tomographic Reconstruction with Tensor Dictionary . . . . . . .
RecT_demo.m . . . . . . . . . . . . . . . . . . . . . . . . .
RecM_Algorithm.m . . . . . . . . . . . . . . . . . . . . . .
2.5.1
2.5.2
2.6.1
2.6.2
3 Alphabetic List of Package Functions
13
13
14
14
14
16
16
16
19
Introduction
1
1.1 Overall
In our recent works [5, 6] we have developed a two-stage algorithm using a dictio-
nary learning approach in a discrete tomographic reconstruction problem both
with a matrix and a tensor formulation. We first build a training data base from
given training images and then construct a dictionary purely from the training
data by looking for a number of basis elements in terms of the dictionary and the
coefficients/representations of the training data in the dictionary. In these two
works ([5, 6]) dictionary learning is formulated as a non-negative matrix/tensor
formulation with sparsity constraints on the representation, and then the recon-
struction problem is formulated as a sparse approximation problem in a convex
optimization framework. We optimized the dictionary learning problem locally
by using the Alternating Direction Method of Multipliers (ADMM), see e.g., [2].
The tomographic reconstruction problems are solved using the software package
TFOCS (Templates for First-Order Conic Solvers) [1].
This documentation describes the implementations details of the proposed
algorithms for solving the tomographic image reconstruction problem in Mat-
lab. This package is designed for a discretized computed tomography problem
formulation, i.e., we use a matrix formulation of the reconstruction problem
(Ax ≈ b), where b ∈ Rm contains the noisy measurement data and A is the
system matrix. The vector x ∈ Rn represents an M × N image, with n = M N,
of absorption coefficients. The total number of tomographic measurement is
m. The system matrix A ∈ Rm×n is available. We use the AIR Tools [4] to
compute the system matrix A. The discretized tomographic problem is a large
sparse and often ill-posed system and incorporating a priori information about
the solution is a necessity to improve the reconstruction and regularize the solu-
tion. This package aims at providing a computational framework for the use of
training images as priors for the solution in tomographic image reconstruction
in the scheme described above. The algorithms work for any size of the sys-
tem matrix A; however, we are more concerned with underdetermined problems
where m < n, because the need for regularization is even more pronounced in
that scenarios.
This package of Matlab routines provides the user with easy-to-use routines
and demo scripts, based on formulation and numerical algorithms described in
[5, 6]. With the demo script files, the user can easily specifies a few values of
problem parameters to obtain a solution to the dictionary learning problem,
2
S. Soltani
A Matlab Package for DL Approach to CT
the tomographic reconstruction problem and compute the best approximation
error of a given test image by the obtained dictionary. For an advanced use
of the routines, we invite the reader to carefully read [5, 6], as well as this
documentation and the descriptions in the .m files provided in this package.
1.2 Toolbox and Software Dependencies
The following toolboxes and software routines are required for this package to
be run smoothly on a computer system.
• Matlab : The package has been implemented in Matlab version 2014a.
• AIR Tools [4] is a Matlab software package for tomographic reconstruc-
tion consisting of tomographic test problems and a number of algebraic
iterative reconstruction methods. The user can download AIR Tools in a
zip format file from http://www2.compute.dtu.dk/~pcha/AIRtools/.
• TFOCS [1] is a suite of programs and routines designed to help construct-
ing of first-order methods for a variety of convex optimization problems
arise in compressed sensing, sparse recovery, and low rank matrix com-
pletion e.g. (BPDN and Lasso). The user can download TFOCS in a zip
format file from http://cvxr.com/tfocs/download/.
The user may download and unpack the software routines in a separate folder
while the paths to those folders can be added to the Matlab’s current path by
including addpath statements.
1.3 Test Images
A collection of gray scale test images has been provided in the package, lo-
cated in the folder “/TestImages”. There are 6 built-in test problems in this
version namely: “Peppers, Matches, Binary, D53, Zirconium and Steel”.
The high-resolution photos of Peppers and Matches are provided by profes-
sor Samuli Siltanen from University of Helsinki. The Zicronium test image is
courtesy of Dr. Hamidreza Abdolvand from University of Oxford. The steel
microstructure image is from [8] and the D53 test image is chosen from the
normalized brodatz texture database [7]. The binary test image is generated by
means of phantomgallery.m function from AIR Tools. The user should update
Matlab’s path with links to the subdirectories /TestImage and /TestResults.
1.4 Notation
The linear tomographic reconstruction problem is formulated as Ax ≈ b where
the vector x ∈ Rn represents the unknown image of size M × N, the vector
b ∈ Rm is the given noisy data, and the matrix A ∈ Rm×n represents the
forward model. We use the following notation, where B is an arbitrary matrix:
,
3
kBkF =P
ij B2
ij
1/2
kBksum =P
ij |Bij|,
S. Soltani
A Matlab Package for DL Approach to CT
2 Package Details
2.1 Dictionary Learning Problem
DL_demo.m
DicLear_M.m
Lea_Algorithm_D2.m
Demo and example script for the dictionary
learning (DL) problem in the matrix form.
The function routine to learn a matrix dictio-
nary with the given parameters.
Performs the ADMM algorithm in k iterations
for the dictionary learning problem in the Ma-
trix form such that D ∈ Ð2.
Lea_Algorithm_Dinf.m Performs the ADMM algorithm in k iterations
for the dictionary learning problem in the Ma-
trix form such that D ∈ Ð∞.
Finds a point in the intersection of two convex
sets by iteratively projecting onto each of the
convex sets.
This function is a tool to illustrate the dictio-
nary elements as images.
montageplot.m∗
dykstra.m∗
2.1.1 DL_demo.m
• Description: Demo script for the dictionary learning (DL) problem in
the matrix form.
• Limit: The test images are predefined in the codes.
• Dependencies: DicLear_M, Lea_Algorithm_D2, Lea_Algorithm_Dinf,
montageplot, dykstra
• Usage: DL_demo
2.1.2 DicLear_M.m
• Description: The function loads the test images, extract training patches
from those images and computes a dictionary with the given parameters
by finding a local optimum to:
min
1
2 kY − DHk2
where kHksum =P
D,H
ij |Hij|,
F + λkHksum
s.t.
D ∈ Ð, H ∈ Rs×t
+ ,
(1)
+ |kdjk∞ ≤ 1}
and Ð2 ≡ {D ∈ Rp×s
Ð∞ ≡ {D ∈ Rp×s
p}.
The patches of the dictionary are of size P × Q with p = P Q. The code
accepts a range of different values for λ.
This function can also run without any input parameter and with the
default values defined in the code.
+ |kdjk2 ≤ √
• Dependencies: Lea_Algorithm_D2, Lea_Algorithm_Dinf, dykstra
4
S. Soltani
A Matlab Package for DL Approach to CT
• Usage:
[ D, H, Norm1H, Norm1HCol, Density, NrNZero, t_end, output ] = ...
DicLear_M(TImage, patch_size, s, Lambda, DSet, t, Outdisplay);
• Inputs:
– TImage: The name of the test image, the test image options are:
Peppers, Matches, Binary, D53, Zirconium, Steel
– patch_size: The training patch size
– s: Number of dictionary elements
– Lambda: The regularization parameter, The code accepts a range of
different values for λ.
– DSet: The name of the compact and convex set which the dictionary
belongs to, the options are: D_inf and D_2.
– t: Number of training patches
– Outdisplay: The display option on the screen, when set to 0 no
display on the screen, when set to 1 prints progress of the algorithm,
default value = 0.
• Outputs:
– D: The matrix dictionary(ies) of size p by s for each λ where
p = patch_size(1)patch_size(2).
– H: The representation matrix(ces) of size s by t for each λ
– Norm1H: kHksum
– Norm1HCol: kHksum/t
– Density: The density percentage(s) of the matrix(ces) H
– NrNZero: The absolute number of nonzeros in the representation
matrix(ces) H
– t_end: Total time used by the ADMM algorithm to find a solution
for a specific λ
– output: Structured output array of the following fields:
F + λkHksum.
∗ Objective: The objective function value(s) for each iteration of
∗ Residual: The residual value(s) for each iteration of the ADMM
∗ StopCr1, StopCr1, StopCr1, StopCr1, StopCr1: The stopping
the ADMM algorithm: 1/2kY − DHk2
algorithm: kY − DHkF.
criteria values at each iteration of the ADMM algorithm.
• Example:
TImage = ’Peppers’;
patch_size = [5 5];
t = 10000;
s = 100;
Lambda = 1;
DSet = ’D_2’;
Outdisplay = 1;
[D,H] = DicLear_M(TImage, patch_size, s, Lambda, DSet, t, Outdisplay);
5
S. Soltani
A Matlab Package for DL Approach to CT
2.1.3 Lea_Algorithm_D2.m
• Description: This function performs the ADMM algorithm in k itera-
tions for the dictionary learning problem in the Matrix form, for one value
of λ, where the dictionary D belongs to the set Ð2, DSet = ’D_2’;
• Usage:
[ W, H, Rs1, Rs2, Rs3, Rs4, Rs5, fobjec, ResF ] = ...
Lea_Algorithm_D2(X, tol, maxiter, s, Cof_Lambda, U, V, ...
Cof_rho, Outdisplay);
• Input:
– X: Training matrix
– tol: The tolerance for the ADMM convergence
– maxiter: Total number of ADMM iterations
– s: The number of dictionary elements
– Cof_Lambda: the regularization parameter λ in the dictionary learn-
ing problem formulation
– U,V: The initial value for the auxiliary variables in the ADMM algo-
rithm
– Cof_rho: The augmented Lagrangian parameter
– Outdisplay: The display option on the screen, when set to 0 no
display on the screen, when set to 1 prints progress of the algorithm,
default value = 0.
• Output:
– W: The dictionary matrix for the corresponding λ
– H: The representation matrix for the corresponding λ
– fobjec: The objective function value for each iteration of the ADMM
2kY − DHk2
F + λkHksum
algorithm: 1
kY − DHkF
– ResF: The residual value for each iteration of the ADMM algorithm:
– Rs1, Rs2, Rs3, Rs4, Rs5: The stopping criteria values at each itera-
tion of the ADMM algorithm
2.1.4 Lea_Algorithm_Dinf.m
• Description: This function performs the ADMM algorithm in k itera-
tions for the dictionary learning problem in the Matrix form„ for one value
of λ, where the dictionary D belongs to the set Ð∞, DSet = ’D_inf’;
• Usage:
[ W, H, Rs1, Rs2, Rs3, Rs4, Rs5, fobjec, ResF ] = ...
Lea_Algorithm_Dinf(X, tol, maxiter, s, Cof_Lambda, U, V, ...
Cof_rho, Outdisplay);
• Input, Output: See Lea_Algorithm_D2.m
6
S. Soltani
A Matlab Package for DL Approach to CT
2.2 Mean Approximation Error
MAEM_demo.m Demo script file illustrating how to compute
the approximation error of a test problem in
the cone defined by the given dictionary.
Computes the approximation error for the ma-
trix formulation.
MAE_M.m
2.2.1 MAEM_demo.m
• Description: This demo script file illustrates how to compute the ap-
proximation error of a test problem in the cone defined by the dictionary.
This code handle various dictionaries for a range of values for λ
• Dependencies: TFOCS Toolbox, MAE_M.m
• Limit: This codes uses TFOCS optimization solver (http://cvxr.com/
tfocs/). Make sure TFOCS is properly installed on the computer.
• Usage: MAEM_demo
2.2.2 MAE_M.m
• Description: This function computes the approximation error for the
matrix formulation i.e., how well we can represent the exact image in the
cone defined by the dictionary. To evaluate the approximation error, i.e.,
the distance of the exact image xexact to its projection on the cone C =
{Dz|z ∈ Rs+} ⊆ Rp
j to the q approximation
problems for all blocks j = 1, 2, . . . , q in xexact. This code handle various
dictionaries for a range of values for λ.
+, we compute the solutions α?
• Dependencies: TFOCS Toolbox
• Limit: The test images are predefined, the user should modify the test
images in the code for other test problem or image sizes. Make sure
TFOCS is properly installed on the computer.
• Usage:
[ Appr_Error_Mean, Appr_Error, Alph, Norm1Alph, SparsityAlph ] =
MAE_M( TrImage, D, Lambda, patch_size, Outdisplay );
• Inputs:
– TrImage: The test image, options are: ’Peppers’, ’Matches’,
’Binary’, ’Zirconium’, ’Steel’, ’D53’
– D: The given(learned) dictionary, note that the dictionary should be
obtained from a training image similar to the test image
– Lambda: The regularization parameter in the dictionary learning
problem formulation, note that it should be consistence with the
dictionary D
– patch_size: The patch sizes of the dictionary elements, should be
provided if the dictionaries patch sizes are rectangular,
default = [ceil(sqrt(size(D,1))),ceil(sqrt(size(D,1)))]
7
S. Soltani
A Matlab Package for DL Approach to CT
– Outdisplay: The display option on the screen, when set to 0 no
display on the screen, when set to 1 prints progress of the algorithm,
default value = 0.
• Outputs:
– Appr_Error_Mean: Mean approximation error
– Appr_Error: The vector of approximation errors for each block in
the image
– Alph: The best representation/approximation of the jth block in the
– Norm1Alph: kα?
– SparsityAlph: The sparsity percentage of the representation vector
jk1 for each block j
cone (α?
j).
j for each block
α?
• Example:
TImage = ’Peppers’;
load(’Dic.mat’,’D’);
Lambda = 1;
patch_size = [5 5];
Outdisplay = 1;
[ Appr_Error_Mean, Appr_Error ] = ...
MAE_M( TrImage, D, Lambda, patch_size, Outdisplay );
2.3 The Reconstruction Problem
RecM_demo
RecM_Algorithm
Perm_Vec∗
L_Matrix∗
Linear_Opr∗
Demo script for the tomographic reconstruc-
tion problem in the matrix form.
Solves the tomographic reconstruction prob-
lem in the matrix form for the asked test im-
age by TFOCS.
Returns the permutation vector which gives a
permutation to the solution.
Returns the matrix L, used to penalize the
block artifacts in the reconstruction formula-
tion.
Defines a linear operator for the tomographic
reconstruction matrix formulation used by the
TFOCS.
2.3.1 RecM_demo.m
• Description: Demo script for the tomographic reconstruction problem in
the matrix form. This script illustrates the use of the dictionary learning
approach in the discrete tomographic reconstruction problem. Note that
this script solves a large scale sparse approximation problem and many it-
erations are needed to converge to the solution and this is not an bug/error
of the code.
• Dependencies: TFOCS, RecM_Algorithm, Perm_Vec, L_Matrix,
Linear_Opr
8