Algorithms for Polyenergetic Breast Tomosynthesis Image - - PowerPoint PPT Presentation

Motivation Mathematical Formulation Reconstruction and Results Final Remarks

Algorithms for Polyenergetic Breast Tomosynthesis Image Reconstruction

Julianne Chung∗, James Nagy†, Ioannis Sechopoulos‡

∗Department of Computer Science

University of Maryland, College Park, MD

†Department of Mathematics and Computer Science

Emory University, Atlanta, GA

‡Department of Radiology and Winship Cancer Institute

Emory University, Atlanta, GA

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks

Outline

1

Motivation

2

Mathematical Formulation

3

Reconstruction and Results

4

Final Remarks

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Outline

1

Motivation

2

Mathematical Formulation

3

Reconstruction and Results

4

Final Remarks

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Breast Cancer Statistics

2nd most common type

• f cancer worldwide

> 1 million women each year diagnosed with breast cancer ∼ 465,000 die each year from the disease

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Digital Tomosynthesis

X-ray Mammography Digital Tomosynthesis Computed Tomography

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

An Inverse Problem

Given: 2D projection images Goal: Reconstruct a 3D volume True Images

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Simulated Problem

Original object: 300 × 300 × 200 voxels (7.5 × 7.5 × 5 cm) 21 projection images: 200 × 300 pixels (10 × 15 cm) −30◦ to 30◦, every 3◦ Reconstruction: 150 × 150 × 50 voxels (7.5 × 7.5 × 5 cm)

Detector Center of Rotation Compressed Breast X-ray Tube Support Plate Compression Plate X-ray Tube Chest Wall Detector

Front view Side view with X-ray tube at 0◦

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Polyenergetic Model

Incident X-ray has a distribution

• f energies

Rh/Rh 28kVp: 47 energy levels, 5keV - 28keV Consequences: Beam Hardening: Low energy photons preferentially absorbed Unnecessary radiation Linear attenuation coefficient depends on energy

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Polyenergetic Model

Incident X-ray has a distribution

• f energies

Rh/Rh 28kVp: 47 energy levels, 5keV - 28keV Consequences: Beam Hardening: Low energy photons preferentially absorbed Unnecessary radiation Linear attenuation coefficient depends on energy

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Monoenergetic Algorithm

Lange and Fessler’s Convex MLEM Algorithm Beam hardening artifacts Monoenergetic Reconstruction

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks An Important Application Polyenergetic Model, Monoenergetic Algorithm

Previous Methods

Methods for eliminating beam hardening artifacts: Dual Energy Methods

Alvarez and Macovski (1976), Fessler et al (2002)

FBP + Segmentation

Joseph and Spital (1978)

Filter function based on density

De Man et al (2001), Elbakri and Fessler (2003)

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Outline

1

Motivation

2

Mathematical Formulation

3

Reconstruction and Results

4

Final Remarks

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

A Polyenergetic Mathematical Representation

Energy-dependent Attenuation Coefficient: µ(e)(j) = s(e)x(j) + z(e) Voxel j where x(j) represents unknown glandular fraction of jth voxel s(e) and z(e) are known linear fit coefficients

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Computing Image Projections

Ray Trace:

• Li

µ(e)dl ≈

N

• j=1

µ(e)(j)a(ij) Vector Notation µ(e) = s(e)x+z(e) ⇒ s(e)Aθx+z(e)Aθ1

1 2 3 4 5 6 2 4 6 10 20 30 40 50 60

Li Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Computing Image Projections

Ray Trace:

• Li

µ(e)dl ≈

N

• j=1

µ(e)(j)a(ij) Vector Notation µ(e) = s(e)x+z(e) ⇒ s(e)Aθx+z(e)Aθ1

1 2 3 4 5 6 2 4 6 10 20 30 40 50 60

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Computing Image Projections

Ray Trace:

• Li

µ(e)dl ≈

N

• j=1

µ(e)(j)a(ij) Vector Notation µ(e) = s(e)x+z(e) ⇒ s(e)Aθx+z(e)Aθ1

1 2 3 4 5 6 2 4 6 10 20 30 40 50 60

Polyenergetic Projection:

ne

• e=1

̺(e) exp (−[s(e)Aθxtrue + z(e)Aθ1])

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Statistical Model

Given x, define the expected value, ¯ b(i)

θ , as the ith entry of ne

• e=1

̺(e) exp (−[s(e)Aθx + z(e)Aθ1]) Let ¯ η(i) be the statistical mean of the noise. Then ¯ b(i)

θ + ¯

η(i) ∈ R is the expected or average observation. Observed Data: b(i)

θ

∼ Poisson(¯ b(i)

θ + ¯

η(i))

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Problem Set-up Statistical Assumptions

Statistical Model

Likelihood Function: p(bθ, x) =

M

• i=1

e−(¯

b(i)

θ +¯

η(i))(¯

b(i)

θ + ¯

η(i))b(i)

θ

b(i)

θ !

Negative Log Likelihood Function: −Lθ(x) = − log p(bθ, x) =

M

• i=1

(¯ b(i)

θ + ¯

η(i)) − b(i)

θ log(¯

b(i)

θ + ¯

η(i))

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Outline

1

Motivation

2

Mathematical Formulation

3

Reconstruction and Results

4

Final Remarks

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Volume Reconstruction

Maximum Likelihood Estimate: xMLE = argmin

x

nθ

• θ=1

−Lθ(x)

• Numerical Optimization:

Gradient Descent: xk+1 = xk − αk∇L(xk) Newton Approach: xk+1 = xk − αkH−1

k ∇L(xk)

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Optimization Tools

Gradient: ∇L(xk) = ATvk where

v(i) =

• b(i)

¯ b(i) + ¯ η(i) − 1 ne

• e=1

̺(e)s(e) exp

• −
• s(e)aT

i xk + z(e)aT i 1

• Hessian: Hk = ATWkA

Hksk = −∇L(xk) ⇔ min

sk ||W

1 2

k Ask − W − 1

2

k

v||2

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Numerical Results

Initial guess: 50% glandular tissue Newton-CG inner stopping criteria:

Max 50 inner iterations residual tolerance < 0.1

Gradient Descent Newton-CG Iteration Relative Error Iteration Relative Error CGLS iterations 1.7707 1.7707

• 1

1.0963 1 1.1399 4 5 0.8753 2 0.8502 3 10 0.8322 3 0.8007 7 25 0.8026 4 0.8045 41

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Compare Images

True

20 40 60 80 100

Mono−convex

0.5 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9

Gradient

20 40 60 80 100

Newton−CG

20 40 60 80 100

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks Reconstruction Algorithms Numerical Results Some Considerations

Convexity

Monoenergetic case (w.r.t. µ)

NO noise or scatter ⇒ ML function is convex Noise and scatter ⇒ ML function SEEMS to be unimodal

Polyenergetic case (w.r.t. density)

Severe nonlinearities ⇒ Cost function is not convex

With new polyenergetic formulation, if A is full rank and

b(i)

θ − (¯

b(i)

θ + ¯

η(i)) ≤ mine s(e) maxes(e) b(i)

θ

¯ b(i)

θ + ¯

η(i) ! ¯ b(i)

θ ,

then we have concavity w.r.t. x.

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks

Outline

1

Motivation

2

Mathematical Formulation

3

Reconstruction and Results

4

Final Remarks

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks

Conclusions

Novel mathematical framework for incorporating a polyenergetic spectrum in tomosynthesis Standard optimization made feasible with new derivation Better reconstructed images, reduced beam hardening artifacts Significant application to breast imaging

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction

Motivation Mathematical Formulation Reconstruction and Results Final Remarks

Future Work

Large-Scale Implementations

Precompute ray-trace matrix for specific geometry Efficient preconditioners for the Newton system

Regularization methods and parameters

Need good regularizer:

Huber penalty, Markov Random Fields, Total Variation

Need good methods for choosing λ

Robustness upon encountering materials not in the model Extension to other nonlinear tomographic applications

(e.g. atmospheric emission tomography, ground penetrating radar, ...)

Thank you!

Julianne Chung, James Nagy, Ioannis Sechopoulos Polyenergetic Tomosynthesis Reconstruction