Upgrade to Pro — share decks privately, control downloads, hide ads and more …

Sparse and Deep Learning Approaches for Biomedical Image Reconstruction

Sparse and Deep Learning Approaches for Biomedical Image Reconstruction

Summer School Lecture by Jong Chul Ye, 13th IEEE EMBS International Summer School on Biomedical Imaging, June 14-21, 2018 , Saint-Jacut de la Mer, France

A3d61bc22cd700a92e7d4136a4d29e8f?s=128

Jong Chul Ye

June 14, 2018
Tweet

Transcript

  1. Jong Chul Ye Bio-Imaging & Signal Processing Lab. Dept. Bio

    & Brain Engineering Dept. Mathematical Sciences KAIST, Korea IEEE EMBS Summer School on Biomedical Imaging Sparse and Deep Learning Approaches for Biomedical Image Reconstruction: Part I: Compressed Sensing & Sparse Recovery This material can be downloaded from http://bispl.weebly.com
  2. Bio Imaging Signal Processing Laboratory Prof. Jong Chul Ye jong.ye@kaist.ac.kr

    h#p://bisp.kaist.ac.kr *4.3.SFDPO$IBMMFOHFXJOOFS LU'0$644  /PBMJBTJOH "DDFMFSBUFE 'BTUJNBHJOHVTJOH$PNQSFTTFETFOTJO H .3*"DDFMFSBUJPO /FVSP*NBHJOH TUBUJTUJDBMUPPMEFWFMPQNFOU  .BDIJOF-FBSOJOH%FFQ-FBSOJOHGPS#JPNFEJDJOF /*3441. .3*4JHOBM1SPDFTTJOH .3BSUJGBDUSFNPWBMVTJOH"-0)" 0QUJDBM*NBHJOH 4VQFSSFTPMVUJPOJNBHJOH '"-$0/  %JGGVTFPQUJDBMUPNPHSBQIZ %05  YSBZ$53FDPOTUSVDUJPO 1&53FDPOTUSVDUJPO $PSSFDUFE 4DBUUFSFE %ZOBNJD%1&5SFDPOTUSVDUJPO The ima ge can not be The ima ge can not be The ima ge can not be The ima ge can not be (SPVOEUSVUI TPVSDF PCKFDU '07 EFUFDUPS #BDLQSPKFDUJPO SBOHF $IPSEMJOF Z Y *OUFSJPS5PNPHSBQIZ Bio-signal processing & Machine learning Applied mathematics, signal processing & machine learning for bio-medical imaging ""1.-PX%PTF$5(SBOE$IBMMFOHF8JOOFS 8BWF/FU  4QBSTF41.
  3. Course Overview 1.  Introduction 2.  Part I: Compressive Sensing • 

    Motion-compensated CS •  Learning-based CS •  Multiple measurement vector (MMV) CS 3.  Part II: Low-rank Hankel Matrix Approaches 4.  Part III: Deep learning approaches 5.  Open problems and outlook
  4. Compressed Sensing Hankel Structured Matrix Comple7on Deep Learning 4 Roadmap:

    From CS to Deep Learning Coherent Theme of Sparse Recovery
  5. Resolution Limits in Medical Imaging •  Diffraction Limit –  Resolution

    limit for optics, x-ray, PET, and etc. http://www.microscopyu.com/
  6. Resolution Limits in Medical Imaging •  Nyquist Sampling limit – 

    Nyquist requires the minimun two times sample rate w.r.t. signal bandlimit –  Resolution limit for Fourier based imaging (MRI, etc.)
  7. Law of Parsimony: Sparsity, Low-rank, LDM •  Occam’s Razor: law

    of parsimony by William Occam (14th century) Entities should not be multiplied beyond necessity. •  All things being equal, the simplest solution tends to be the best one.
  8. Super-resolution Microscopy S. Hell et al, Science 2007.

  9. Compressed Sensing MR Imaging •  Forward problem •  Sparse recovery

    problem Figure courtesy of Mathews Jacob CS
  10. 9 View CT Reconstruction using Deep Learning

  11. COMPRESSED SENSING

  12. Compressed Sensing (CS) 12 measurements vector # non-zeros •  Incoherent

    projection •  Underdetermined system •  Sparse unknown vector Courtesy of Dr. Dror Baron b A x n ⇥ 1 k m ⇥ 1 m ' k log( n ) ⌧ n
  13. Min-norm Solution is Not sparse

  14. L0 Minimization •  Two fundamental questions –  Uniqueness ? – 

    Convex relaxation ? min x k x k0 subject to b = A x
  15. Uniqueness and Spark

  16. L1 Convex Relaxation min x k x k1 subject to

    b = A x k x k1 = n X i=1 | xi |
  17. Intuition behind l1 R. Tibshirani, “Regression Shrinkage and Selection via

    the Lasso”, Jour. Royal Stat. Soc. B, 1996,
  18. Sufficient Conditions for P1

  19. RIP Conditions •  RIP has been developed for probabilistic argument

    .
  20. RIP as Sufficient Condition

  21. RIP: Robust Recovery

  22. 22 Sparse Recovery Formulations min x kT x k1 subject

    to k b A x k  ✏ min ✓ k✓k1 subject to kb A ✓k  ✏ •  Analysis Formulation •  Synthesis Formulation
  23. K. Scheffler et al, Eur Radiol (2003) 13:2409–2418 Sequence diagram

    of b-SSFP §  Demanding Acquisition Requirements An example of Cartesian TrueFISP parameter for cardiac imaging. The acquisition time is calculated based on 60bpm patients with 8 slice acquisition. Flip angle of 40-60° is commonly used. §  Balanced SSFP FLASH b-SSFP 23 / 38 Application : Cardiac MRI
  24. Aliasing in accelerated Dynamic MRI t Down sampling along k-t

    space à Aliasing artifact
  25. Classical Methods •  Parallel imaging (SENSE, GRAPPA) –  Prussmann et

    al, 1999, Griswold et al 2002 •  k-t space method –  UNFOLD –  UNFOLD-SMASH, TSENSE –  K-t BLAST/SENSE –  PARADISE, PARADIGM Blaimer et al, Top Magn Reson Imaging • V15, N 4, August 2004
  26. How to Improve ? Sparsity !! t y KLT Cardiac

    MR fMRI
  27. RIP in Dynamic MR k-t space sample Random phase encoding

    Temporal transform k t
  28. k-t FOCUSS (Jung et al, PMB:2007, MRM:2009(a), MRM:2009(b)) k-t FOCUSS

    (our method) k-t BLAST/SENSE (J. Tsao, MRM, 2003) No-update Optimal from compressed sensing perspective Does not solve L1 minimization problem Not Optimal from CS perspective When there is no update, k-t FOCUSS is exactly same with k-t BLAST/SENSE
  29. Zero-padding inverse FFT from measurements 11 x accelerated Sampling pattern

  30. k-t BLAST 11 x accelerated Sampling pattern

  31. k-t FOCUSS with temporal average 11 x accelerated Sampling pattern

  32. MOTION-COMPENSATED COMPRESSED SENSING

  33. Lessons from HDTV History MUSE MPEG The first HDTV (Japan)

    (obsolete) Modern Standard (MPEG1,2,4..) quincunx sampling ME/MC Lattice sampling theory à Low compression ME/MC + residual coding à High compression * MUSE: MUltiple Sub-nyquist Encoding systems Figures from Kovacevic et al, IEEE TIP, 1993
  34. MUSE vs k-t MR MUSE k-t MR The first HDTV

    (Japan) (obsolete) UNFOLD k-t BLAST/SENSE quincunx sampling k-t space downsampling Lattice sampling theory à Low compression Lattice sampling theory à Low acceleration !! k time * MUSE: MUltiple Sub-nyquist Encoding systems Figures from Kovacevic et al, IEEE TIP, 1993
  35. Reference is often Free •  Cardiac MR –  During volume

    acquisition, diastole phase can be acquired •  fMRI •  MR angiography, vocal tract, etc. •  example) RIGR, SPEAR, etc. 12 sec X 10 Finger Tapping 30 sec 30 sec 21 sec Rest Rest Dummy
  36. Random Sampling with Reference When two reference frames can be

    obtained, k Reference frames Fully sampled center region t Random sampling When one reference frames can be obtained 0 x k t We want to distribute k-space samples to allow prediction and residual coding
  37. ME/MC in MPEG Video t k Encoder Full sampled measurements

    t x y
  38. ME/MC in MPEG Video t x y

  39. ME/MC in MPEG Video t x y Residual encoding

  40. ME/MC in Dynamic MRI t k Encoder Down sampled measurements

    t x y
  41. Refinable Motion Estimation t k Encoder Down sampled measurements t

    x y
  42. Refinable Motion Estimation t x y Initialization with k-t FOCUSS

    k t
  43. Refinable Motion Estimation t x y Initialization with k-t FOCUSS

    k t
  44. Refinable Motion Estimation t x y Apply ME/MC k t

  45. Refinable Motion Estimation t x y Residual Encoding using k-t

    FOCUSS Residual Encoding ME/MC
  46. Cine MRI •  Acquisition parameters - 1.5 T Philips scanner

    at Yonsei Unversity Medical Center in Korea - bSSFP - FOV = 345.00 x 270.00 mm2 - Matrix size : 256 x 220 ( Read-Out x PE ) - TR = 3.45ms - flip angle : 50° - heart rate : 66 bpm - # of cardiac phases : 25 frames •  Among fully sampled k-data, partial k-data were arbitrarily chosen.
  47. Zero-padding inverse FFT from measurements 11 x accelerated Sampling pattern

  48. k-t BLAST 11 x accelerated Sampling pattern

  49. k-t FOCUSS 11 x accelerated Sampling pattern

  50. Motion Compensated k-t FOCUSS 11 x accelerated Sampling pattern

  51. Radial Acquisition

  52. World-first In vivo Experiment (Joint work with K. Nayak at

    USC) •  Acquisition parameters - 3.0 T GE Signa EXCITE at the University of Southern California - bSSFP - FOV = 300 x 300 mm2 - Matrix size : 240 x 256 ( PE x Read-Out) - TR = 3.7ms -  views per segment = 10 - flip angle : 45° - heart rate : 65 bpm - # of ardiac phases : 25 frames -  acquisition time for 6 x accelerated data = 4 heart beats -  acquisition time for fully sampled data = 24 heart beats
  53. In Vivo Experiment (6x accel.)

  54. LEARNING-BASED COMPRESSED SENSING

  55. Which One Is Optimal ? FT along t-axis PCA along

    t-axis WVT along x-axis + FT along t-axis x t x f
  56. KLT (PCA): optimal transform Karhunen Loeve Transform is well known

    for an optimal sparsifying transform of general signals. KLT basis: Matrix of Eigen vectors for T X U V = Σ 0 0 L s s ⎛ ⎞ ⎜ ⎟ Σ = ⎜ ⎟ ⎜ ⎟ ⎝ ⎠ T T U X V = Σ Considering SVD L: Large Singular values s: small Singular values Hence, KLT is SPARSE 0 L s >> → T XX X
  57. PCA in k-t FOCUSS k-t space sample Random phase encoding

    Temporal transform t y SVD Significant eigenvector
  58. 58 Two step Boosting for PCA basis

  59. MR Angiography Golden Section Radial Trajetory kx ky kx ky

    kx ky kx ky … T=1 T=2 T=512 … ky kx ky kx ky kx … … Target time points: 39, 79, 118, 157, 197, 236, 275, 314, 354, 393, 432, 472 T=39 T=472 39 or 40 views are merged on each time frame Measurement size (512 read-out x 512 views) vs Total recon matrix (512 x 512 x 12 time frames) 12 x acceleration !!
  60. KLT for CE-MRA x t k-t FOCUSS recon using FFT

    Summation along t-axis Thresholding
  61. 2009 ISMRM Recon Challenge Award L2 minimization k-t FOCUSS using

    KLT
  62. 62 T2 Molecular motion Fast (water) Slow (protein) (fat) Detect

    T2 changes due to variety of disease T2 Parameter Mapping (joint work D. Kim, NYU: Feng et al, MRM, 2011)
  63. 63 Cardiac T2 Mapping 6 x accel. Conventional method 6

    x accel. k-t FOCUSS
  64. 64 Cardiac T2 Mapping 1.8 x accel. GRAPPA 6 x

    accel. k-t FOCUSS
  65. 65 T2 Estimation

  66. MULTIPLE MEASUREMENT VECTOR COMPRESSED SENSING

  67. MMV problem - definition minimize kXk 0 subject to B

    = AX. k = kXk 0 : number of nonzero rows r : number of snapshots m : number of sensor elements Here, we assume that the columns of B are linearly independent. Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 3 / 21 MMV: Multiple Measurement Vector Problems
  68. Geophysical x x x x x x x x x

    x x x x x x x x x x Medical x x x x x x Industrial Electromagnetic Acoustic Ultrasonic Optical x x x x x x Electromagnetic Ultrasonic Optical Off-set VSP/ cross-well tomography GPR surface imaging induction imaging Ultrasound tomography optical microscopy photon imaging Ultrasound tomography optical microscopy induction imaging Inverse Scattering Problems Slide courtesy by Devaney
  69. Diffuse Optical Tomography Near-infrared (NIR, ~650-950 nm) * Durduran et

    al., MICCAI, 2010 ⌦1 ⌦2
  70. Elastic Source Imaging Elastography: Medical, geophysical applications https://iame.com/online/breast_elastography/ https://marcellusdrilling.com/2018/04/9-more-seismic- testing-devices-stolen-in-swpa-6-were-returned/

  71. Joint Sparsity Multiple Measurement Vectors from Multiple illumination Patterns Finite

    # of snapshots are available
  72. 72 Sparse Signal Perturbations in the optical parameters Angiogenesis in

    cancer Joint Sparsity
  73. 73 Inverse Scattering Problems Ill-posed Non-linear Existing methods •  Born

    approximation à linearization error •  Iterative Born approximation à computationally expensive
  74. Overcoming Nonlinearities using Joint Sparsity 74

  75. l0 uniqueness results on MMV Definition Given a matrix A,

    let spark(A) denote the smallest number of linearly dependent columns of A. Theorem (Rao / Chen, Huo / Feng, Bresler / Davies, Eldar) If a matrix X satisfies AX = B, then kXk 0 < spark(A) + rank(B) 1 2 if and only if X is the unique solution to the l 0 minimization problem. With increasing number of snapshots, more non-zero elements can be recovered. Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 5 / 21 L0 Uniqueness of MMV
  76. Spark reduction 1 J. Kim, O. Lee, J. Ye, Compressive

    MUSIC (CS-MUSIC), IEEE Trans. on Inform. Theory, January 2012 Theorem Let r  m < n. Suppose that we are given a sensing matrix A 2 Rm⇥n and an observation matrix B 2 Rm⇥r . If the k nonzero rows of X are in general position (i.e., any collection of r nonzero rows are linearly independent) and A satisfies the RIP condition 0  L 2 k r +1 (A) < 1, then spark(Q⇤A) = k r + 1. Note that A 2 Rm⇥n satisfies RIP with 0  L 2 k r +1 (A) < 1 if and only if k < spark(A) + rank(B) 1 2 . :l 0 bound for MMV problem CS-MUSIC achieves l 0 bound with finite snapshot. Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 6 / 21 Spark Reduction Principle
  77. Generalized MUSIC criterion Theorem Assume that we have an MMV

    problem AX = B, where A, X and B as in the previous theorem. If Ik r ⇢ suppX with |Ik r | = k r and AIk r 2 Rm⇥ ( k r ) which consists of columns of A, whose indices are in Ik r , then for any j 2 {1, · · · , n} \ Ik r , a⇤ j h PR ( Q ) PR ( PR ( Q ) AIk r ) i aj = 0 () j 2 suppX. a⇤ j P? [ AIk r B ] aj = 0 () j 2 suppX. Augmented Signal Subspace: R[AIk r B] Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 7 / 21 Generalized MUSIC Criterion
  78. Geometric Interpretation Jong Chul Ye (KAIST) The role of joint

    sparsity UIUC Talk 8 / 21 Geometric Interpretation
  79. Su cient condition for SS-OMP to find a partial support

    Case 1: r is a fixed number Theorem Let k ⇣ [ AIk r S ] ⌘ > 1 + k r . If we have m > 1 + 1 2k r ⌘k r 2k log (n k) r , then we can find k r correct indices of suppX by applying subspace S-OMP. 1. [Fletcher, Rangan] For SMV, when SNR ! 1, we need m > 2(1 + )k log (n k) for some > 0. 2. The sampling ratio is reciprocally proportional with respect to the number of multiple measurement vectors. 3. SNR! 1 is required, in the large system limit. Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 9 / 21 (Similar to Reevs and Gaspar) Sampling Rate for MMV
  80. Su cient condition for S-OMP to find a partial support

    Case 2: ↵ := limn!1 r(n)/k(n) > 0 Theorem Let k ⇣ [ AIk r S ] ⌘ > 1 + ↵ Then if we have m > (1 + )2 k 1 2⌘k r ↵ [2 F(↵)]2 , for some > 0 where F(↵) is an increasing function on [0, 1] such that F(1) = 1 and lim↵! 0+ F(↵) = 0. Then we can find k r correct indices of suppX by applying subspace S-OMP. 1. If ↵ ! 1 and SNR ! 1, then we only need to have m > (1 + )k for a small > 0. (cf. MUSIC) 2. In this case, the required number of sensor elements is at most 4(1 + )k when SNR! 1 (No log n factor). 3. The required SNR is finite. Jong Chul Ye (KAIST) The role of joint sparsity UIUC Talk 10 / 21 Sampling Rate for MMV
  81. Other MMV Applications •  Parallel MRI + CS •  EEG/MEG

    •  Diffuse optical tomography •  Wave inverse scattering MMV for Bio-medical Imaging IFFT MR sca n
  82. Optical Diffraction Tomography (ODT) Laser L1 L2 P BS1 GM

    L3 CL OL L4 L5 M1 L6 M2 BS2 Camera M3 P
  83. ODT under Born/Rytov Approximation Sample z x(y) 3d k-space Scattered

    field 2 2-D Fourier transform and scaling
  84. ODT using Joint Sparsity Lippman-Schwinger equation MMV for Support &

    Induced Current Recovery
  85. ODT using Joint Sparsity Unknown Total Field Estimation Abnormality Estimation

  86. Super-Resolution from Joint Sparse Recovery (Lim et al, OPEX, 2017)

  87. Super-Resolution from Joint Sparse Recovery (Lim et al, OPEX, 2017)

      = ; < ; < = = ; < ; < = D 5\WRY E 3URSRVHG
  88. Super-Resolution from Joint Sparse Recovery (Lim et al, OPEX, 2017)

    = ; < ; < = = ; < ; < =   D 5\WRY E 3URSRVHG
  89. 89 Diffuse Optical Tomography Applications (Lee, OPEX, 2013)

  90. 90 Diffuse Optical Tomography Applications (Lee, OPEX, 2013)

  91. Summary : Part I •  Sparsity principle is very important

    for biomedical image reconstruction︎ •  Compressed sensing has many extensions︎ •  Compressed sensing︎ •  Motion compensated CS︎ •  Learning-based CS︎ •  MMV CS︎ •  Downside: suffers from discretization, RIP, optimization bottleneck︎
  92. References: Part I •  Ye, J. C., Tak, S., Han,

    Y., & Park, H. W. (2007). Projection reconstruction MR imaging using FOCUSS. Magnetic Resonance in Medicine, 57(4), 764-775. •  Jung, H., Ye, J. C., & Kim, E. Y. (2007). Improved k–t BLAST and k–t SENSE using FOCUSS.Physics in Medicine and Biology, 52(11), 320 •  Ye, J. C. (2007). Compressed sensing shape estimation of star-shaped objects in Fourier imaging. IEEE Signal Processing Letters, 14(10), 750-753. •  Jung, H., Sung, K., Nayak, K. S., Kim, E. Y., & Ye, J. C. (2009). K-t FOCUSS: A general compressed sensing framework for high resolution dynamic MRI. Magnetic Resonance in Medicine, 61(1), 103-116 •  Jin, K. H., Kim, Y., Yee, D. S., Lee, O. K., & Ye, J. C. (2009). Compressed sensing pulse-echo mode terahertz reflectance tomography. Optics Letters, 34(24), 3863-3865. •  Jung, H., & Ye, J. C. (2010). Motion estimated and compensated compressed sensing dynamic magnetic resonance imaging: What we can learn from video compression techniques. International Journal of Imaging Systems and Technology, 20(2), 81-98. •  Jung, H., Park, J., Yoo, J., & Ye, J. C. (2010). Radial k-t FOCUSS for high‐ resolution cardiac cine MRI. Magnetic Resonance in Medicine, 63(1), 68-7
  93. References: Part I •  Choi, J., Kim, K. S., Kim,

    M. W., Seong, W., & Ye, J. C. (2011). Sparsity driven metal part reconstruction for artifact removal in dental CT. Journal of X-ray Science and Technology, 19(4), 457-475. •  Lee, O., Kim, J. M., Bresler, Y., & Ye, J. C. (2011). Compressive diffuse optical tomography: noniterative exact reconstruction using joint sparsity. IEEE Transactions on Medical Imaging, 30(5), 1129-1142. •  Lee, K., Tak, S., & Ye, J. C. (2011). A data-driven sparse GLM for fMRI analysis using sparse dictionary learning with MDL criterion. IEEE Transactions on Medical Imaging, 30(5), 1076-1089. •  Feng, L., Otazo, R., Jung, H., Jensen, J. H., Ye, J. C., Sodickson, D. K., & Kim, D. (2011). Accelerated cardiac T2 mapping using breath-hold multiecho fast spin- echo pulse sequence with kt FOCUSS. Magnetic Resonance in Medicine, 65(6), 1661-1669. •  Kim, J. M., Lee, O. K., & Ye, J. C. (2012). Compressive MUSIC: Revisiting the link between compressive sensing and array signal processing. IEEE Transactions on Information Theory, 58(1), 278-301.
  94. References: Part I •  Kim, J. M., Lee, O. K.,

    & Ye, J. C. (2012). Improving noise robustness in subspace-based joint sparse recovery. IEEE Transactions on Signal Processing, 60(11), 5799-5809. •  Min, J., Jang, J., Keum, D., Ryu, S. W., Choi, C., Jeong, K. H., & Ye, J. C. (2013). Fluorescent microscopy beyond diffraction limits using speckle illumination and joint support recovery. Scientific Reports, 3(2075) •  Lee, O., & Ye, J. C. (2013). Joint sparsity-driven non-iterative simultaneous reconstruction of absorption and scattering in diffuse optical tomography. Optics Express, 21(22), 26589-26604. •  Zong, X., Lee, J., Poplawsky, A. J., Kim, S. G., & Ye, J. C. (2014). Compressed sensing fMRI using gradient-recalled echo and EPI sequences. NeuroImage, 92, 312-321. •  J. Min C. Vonesch, H.Kirshner, L. Carlini, N. Olivier, S. Holden, S. Manley, J.C. Ye, M. Unser, (2014) FALCON: fast and unbiased reconstruction of high- density super-resolution microscopy data, Scientific Reports 4 , Article no 4577.
  95. References: Part I •  Yoon, H., Kim, K. S., Kim,

    D., Bresler, Y., & Ye, J. C. (2014). Motion adaptive patch-based lowrank approach for compressed sensing cardiac cine MRI. IEEE Transactions on Medical Imaging, 33(11), 2069-2085. •  Min, J., Holden, S. J., Carlini, L., Unser, M., Manley, S., & Ye, J. C. (2014). 3D high-density localization microscopy using hybrid astigmatic/biplane imaging and sparse image reconstruction. Biomedical Optics Express, 5(11), 3935-3948. •  Han, P. K., Park, S. H., Kim, S. G., & Ye, J. C. (2015). Compressed sensing for fMRI: feasibility study on the acceleration of non-EPI fMRI at 9.4 T. BioMed research international, 2015. •  Ward, J. P., Lee, M., Ye, J. C., & Unser, M. (2015). Interior tomography using 1D generalized total variation. Part I: Mathematical foundation. SIAM Journal on Imaging Sciences, 8(1), 226-247 •  Lee, M., Han, Y., Ward, J. P., Unser, M., & Ye, J. C. (2015). Interior tomography using 1D generalized total variation. Part II: Multiscale implementation. SIAM Journal on Imaging Sciences, 8(4), 2452-2486.
  96. References: Part I •  Kim, K., Ye, J. C., Worstell,

    W., Ouyang, J., Rakvongthai, Y., El Fakhri, G., & Li, Q. (2015). Sparse-view spectral CT reconstruction using spectral patch-based low-rank penalty. IEEE Transactions on Medical Imaging, 34(3), 748-760.. •  Kim, K., Son, Y. D., Bresler, Y., Cho, Z. H., Ra, J. B., & Ye, J. C. (2015). Dynamic PET reconstruction using temporal patch-based low rank penalty for ROI-based brain kinetic analysis. Physics in Medicine and Biology, 60(5), 2019. •  Lee, O. K., Kang, H., Ye, J. C., & Lim, M. (2015). A non-iterative method for the electrical impedance tomography based on joint sparse recovery. Inverse Problems, 31(7), 075002. •  Lim, J., Lee, K., Jin, K. H., Shin, S., Lee, S., Park, Y., & Ye, J. C. (2015). Comparative study of iterative reconstruction algorithms for missing cone problems in optical diffraction tomography. Optics Express, 23(13), 16933-16948. •  Lee, O., Tak, S., & Ye, J. C. (2015). A unified sparse recovery and inference framework for functional diffuse optical tomography using random effect model. IEEE Transactions on Medical Imaging, 34(7), 1602-1615.
  97. References: Part I •  Ye, J. C., Kim, J. M.,

    & Bresler, Y. (2015). Improving M-SBL for joint sparse recovery using a subspace penalty. IEEE Transactions on Signal Processing, 63(24), 6595-6605. •  Young-Beom Lee, Jeonghyeon Lee, Sungho Tak, Kangjoo Lee, Duk L. Na, Sangwon Seo, Yong Jeong, and Jong Chul Ye (2016). Sparse SPM: Group sparse-dictionary learning in SPM framework for resting-state functional connectivity MRI analysis, NeuroImage, 125, 1032–1045 •  Yoo, J., Jung, Y., Lim, M., Ye, J. C., & Wahab, A. (2017). A joint sparse recovery framework for accurate reconstruction of inclusions in elastic media. SIAM Journal on Imaging Sciences, 10(3), 1104-1138.
  98. Jong Chul Ye Bio-Imaging & Signal Processing Lab. Dept. Bio

    & Brain Engineering Dept. Mathematical Sciences KAIST, Korea IEEE EMBS Summer School on Biomedical Imaging Sparse and Deep Learning Approaches for Biomedical Image Reconstruction: Part 2: Low-Rank Hankel Matrix Approaches
  99. Compressed Sensing Hankel Structured Matrix Comple7on Deep Learning 99 Roadmap:

    From CS to Deep Learning Coherent Theme of Sparse Recovery
  100. THEORY

  101. Compressed Sensing (CS) 101 measurements vector # non-zeros •  Incoherent

    projection •  Underdetermined system •  Sparse unknown vector Courtesy of Dr. Dror Baron b A x n ⇥ 1 k m ⇥ 1 m ' k log( n ) ⌧ n
  102. Limitations of CS •  Discrete domain theory •  Explicit form

    of sensing matrix •  RIP issue
  103. Analytic Reconstruction (b) Time-reversal of a scattered wave (a) MR

    Imaging Beautiful analytic reconstruction results from fully sampled data
  104. None
  105. Finite Rate of Innovation Sampling Theory (Vetterlie, 2002) T -T

    0 n1 -n1 0
  106. FRI to Low-Rank Hankel matrix ︙ ︙ 1 2 3

    4 5 6 7 8 9 -1 0 1 2 3 4 5 6 7 8 9 1 2 3 4 5 6 7 8 9 12 2 3 4 5 6 7 8 9 12 13 3 4 5 6 7 8 9 10 10 10 10 0 11 11 11 1 2 3 4 5 Finite length convolution Matrix Representation * ALOHA : Annihilating filter based LOw rank Hankel matrix Approach * Jin KH et al. IEEE TCI,2016 * Jin KH et al.,IEEE TIP, 2015 * Ye JC et al. IEEE TIT, 2016
  107. Hankel Matrix

  108. Low-Rank Hankel matrix minimization Missing elements can be found by

    low rank Hankel structured matrix compleAon Nuclear norm Projec;on on sampling posi;ons min m kH ( m ) k⇤ subject to P⌦(b) = P⌦( f ) RankH(f) = k * Jin KH et al IEEE TCI, 2016 * Jin KH et al.,IEEE TIP, 2015 * Ye JC et al., IEEE TIT, 2016 m
  109. General TV Signals 109 Weighted Fourier data Piecewise smooth Splines,

    polynomials
  110. 110 Existence of Annihilating Filter ˆ h(!) ⇤ ⇣ ˆ

    l(!) ˆ f(!) ⌘ = 0 Annihilating filter for weighted Fourier data General Low-Rank Hankel Matrix Completion
  111. Visualization of concept

  112. Extension to general signal models Stream of Diracs Stream of

    differentiated Diracs Non-uniform spline Piecewise smooth polynomial rank  X j dj rank = r rank  rq rank = r rank = X j dj rank = rq * Ye JC et al.,IEEE TIT 2016 With a proper weighAng, the Hankel matrix of the weighted k-space data à low ranked.
  113. None
  114. Performance Guarantees m c1µcsk log ↵ n min m kH(m)k⇤

    subject to P⌦( m ) = P⌦(ˆ f ) min m kH(m)k⇤ subject to kP⌦( m ) P⌦(ˆ f ) k  kH(m) H(ˆ f)kF  c2n2 ↵ = ( 2 , on grid 4 , o↵ grid Exact Recovery Stable Recovery * Ye JC et al.,IEEE TIT 2016
  115. Low-Rank Hankel matrix minimization Missing elements can be found by

    low rank Hankel structured matrix compleAon Nuclear norm Projec;on on sampling posi;ons min m kH ( m ) k⇤ subject to P⌦(b) = P⌦( f ) m
  116. Optimization 116 •  Complexity depends on an estimated rank of

    patch !! ADMM
  117. APPLICATIONS

  118. Applications to Image Processing Inpainting & Impulse noise removal

  119. Spectral Domain Sparsity 119 Smoothness, texture, pattern Sparse spectrum Smooth

    patch Edge patch Texture patch Concentrated at DC Elongated along Perpendicular to edge Distributed orthogonally w.r.t. texture
  120. Low-rank Hankel matrix in Image

  121. 2-D Hankel matrix 121 Hankel Matrix construction

  122. 18. APR. 2015. 122 Why patch processing ? -  Spectrum

    changes for each patch -  Need to adapt the local Image statistics ALOHA
  123. Rotation invariant sparsity 18. APR. 2015. 123 Hankel structured matrix

    is intrinsic low rank !! 0 deg. 90 deg. 0 90
  124. Experimental results (x5) 18. APR. 2015. 124 *

  125. Experimental results (x5) 18. APR. 2015. 125 *

  126. 18. APR. 2015. 126 Text inlayed image reconstruction

  127. 18. APR. 2015. 127 Line scratches

  128. 18. APR. 2015. 128 Object removal

  129. 129 *K.H. Jin, et. al, IEEE TIP (2015)

  130. Lissajous scanning in MEMS lens scanning OCT microscopy 18. APR.

    2015. 130 Park, Hyeon-Cheol, et al. "Forward imaging OCT endoscopic catheter based on MEMS lens scanning." Optics letters 37.13 (2012): 2673-2675. (*courtesy of H. Park and K. Jeong) * Lissajous scanning •  Lissajous scanning utilizes full speeds of both scanning axes! •  Lissajous scanner in MEMS controlled with only small voltage.
  131. Temporal sampling in Lissajous scanning 18. APR. 2015. 131 Subsampled

    image in the regular grid can be interpolated by inpainting algorithms. Subsampled
  132. Retrospective Experiments 18. APR. 2015. 132 * Data acquired by

    Lina Carlini (LEB, EPFL) *Endoplasmic Reticulum labeled with tdEOS fused to reticulon-4 in U2OS cells was taken on a modified epi-fluorescence microscope (Olympus IX71).
  133. Real OCT experiment (x6) 18. APR. 2015. 133 * Data

    acquired by Hyeon-Cheol Park (Johns Hopkins Univ.) and Ki-Hun Jeong (KAIST) Korean coin 1.  140 x 140 px2 2.  Total data length : 150k
  134. Magnified views (x6) 18. APR. 2015. 134

  135. 135 Impulse Noise Removal (Jin et al, TIP, 2018)

  136. 70% Salt and Pepper Noise Removal

  137. 137 Multi-Channel Impulse Noise Removal

  138. Recall: Sparsity in dynamic MRI t y KLT Cardiac MR

    fMRI
  139. Sparsity of Dynamic MRI y (PE) Sparse In Wavelet domain

    x (RO) y (PE) x (RO) Sparsity #1 Sparsity #2 T2 images s-t spectrum y (PE) t y (PE) f Temporal Fourier tra nsform Hankel matrix with wavelet weigthing ky Hankel matrix on t 2D block Hankel matrix on ky -t t ky
  140. k-t dynamic cases

  141. Dynamic MRI – multi coil Six fold (x6) down sampling

    # of coils =4
  142. Visualization of concept

  143. Single coil static MRI

  144. Parallel MRI Signal Model 144 i-th coil measurement Signals from

    parallel MRI i-th coil sensitivities Original signal × =
  145. Multi-channel Generalized ALOHA in kx -ky 145 FT × ×

    × = × Annihilation #2 : inter-coil relationships
  146. Rank Bound for Parallel Imaging Jin et al, IEEE TCI,

    2016 Sparsity of common image In transform domain Sparsity of sensitivity map In Fourier domain
  147.          

                                                                           DIBOOFM TJHOBM                                                                               DIBOOFM TJHOBM DIBOOFM TJHOBM E E E -PX3BOL .BUSJY $PNQMFUJPO -PX3BOL .BUSJY $PNQMFUJPO -PX3BOL .BUSJY $PNQMFUJPO -1 -1 -1                                                                                                                                                            
  148. Parallel MRI

  149. What is MR parameter mapping? MR Parameter Mapping [1] Siemonsen,

    S et al., Radiology, 2008 [2] Rugg-Gunn, F. J., et al. NEUROLOGY, 2005 TE1 TE2 TE3 TE4 e.g. Multi-Echo Spin-Echo (ME-SE, T2 mapping) t Finding the quanAtaAve v alue of each ;ssue Cons Pros e.g. of T2 mapping for diagnosis for epilepsy Clinically valuable - As a quantitative diagnosis tool - Acute stroke, epilepsy, etc. Pros Cons TE1 TE2 TE3 TE4 e.g. Multi-echo images for T2 mapping Long scan time - Needs multiple scans - Variation of TI, TE, FA, etc.
  150. Sparsity of Dynamic MRI y (PE) Sparse In Wavelet domain

    x (RO) y (PE) x (RO) Sparsity #1 Sparsity #2 T2 images s-t spectrum y (PE) t y (PE) f Temporal Fourier tra nsform Hankel matrix with wavelet weigthing ky Hankel matrix on t 2D block Hankel matrix on ky -t t ky
  151. 47.5173×10-3 52.9943×10-3 47.6076×10-3 47.9392×10-3 17.6058×10-3 60.6220×10-3 ME-SE image Error ×5

    Reconstruction of x12.8 accelerated scan – ME-SE (4th echo) Result : in vivo acceleration study ( ME-SE, T2 ) k-t FOCUSS k-t SPARSE C-LORAKS Full k-t SLR ALOHA Patch Lee et al, MRM, 2015
  152. 18.0169×10-2 15.9958×10-2 18.7474×10-2 23.0070×10-2 6.5596×10-2 12.6578×10-2 0 (ms) 50 100

    150 T2 Map Error ×5 Mapping from the reconstruction of x12.8 accelerated scan – T2 mapping Result : in vivo acceleration study ( ME-SE, T2 ) k-t FOCUSS k-t SPARSE C-LORAKS Full k-t SLR ALOHA Patch Lee et al, MRM, 2015
  153. Result : Signal intensity curves ( SE-IR, T1 ) 400

    800 1200 1600 0 Time (ms) -1 0 1 0SJHJOBM LU'0$644 LU41"34& 1BUDI LU4-3 -03",4 "-0)" Signal Intensity (a.u.) T1 relaxation LU'0$644 LU41"34& 1BUDI 0SJHJOBM "-0)" -03",4 LU4-3 Lee et al, MRM, 2015
  154. Summary of Results •  Goal : Acceleration of MR Parameter

    mapping by undersampling and reconstruction Full acquisition Scan Time Conventional scan Accelerated sc an 1min 2s 12min 50s Accelerated acq uisition Reconstructed i mages ALOHA Lee et al, MRM, 2015
  155. MR artifacts •  What is MR artifacts? During acquisition, external

    interruptions (ex. fluctuation power supply of gradient, motion of object, etc.) distort signals. 155 4QJLFOPJTF 3FTQJSBUPSZNPUJPO M. Graves, et. al., JMRI (2013)
  156. Motivation 156

  157. ü  )FSSJOHCPOF TQJLFT%LTQBDF  Motivations: MR artifacts as sparse outliers

    157
  158.  ü  .PUJPOBSUJGBDU TQJLFT%LTQBDFQBSBMMFMUPSFBEPVU   158 Motivations: MR artifacts

    as sparse outliers
  159. ;JQQFSBSUJGBDU TQJLFT%LTQBDFQFSQFOEJDVMBSUPSFBEPVU  159 Motivations: MR artifacts as sparse outliers

  160. Key Observation : Sparse outliers 160 ALOHA* Sparse MR artifact

    + 4QBSTFPVUMJFSJTTUJMMTQBSTFJOXFJHIUFE)BOLFMNBUSJY † E. Candes, et. al, JACM (2011), R. Otazo, et. al, MRM (2015) *K.H. Jin, et. al, IEEE TIP (2015), K.H. Jin, et. al, arXiv (2015), J. C. Ye, et. al, arXiv (2015), J. Lee, et. al, MRM (2016), D. Lee, et. al, MRM (2016) •  ALOHA: Annihilating filter based LOw rank Hankel matrix Approach
  161. Robust ALOHA 161 RPCA for weighted Hankel matrix ü  Extension

    of ALOHA for decomposition of sparse outliers (E) out of mixed signal* ü  Can be addressed ADMM† ü  K-space weighting signal Sparse outlier K-space weighting &$BOEFT FUBM +"$.  30UB[P FUBM .3.   õ4#PZE FUBM 'PVOEBUJPOTBOE5SFOETJO.BDIJOF-FBSOJOH   ö;8BO FUBM .BUIFNBUJDBM1SPHSBNNJOH$PNQVUBUJPO  
  162. 162 Algorithm Flowchart

  163. Retrospective results )JHIJOUFOTJUZ 4QJLFOPJTF -PXJOUFOTJUZ 4QJLFOPJTF MPXGSFRVFODZSFHJPO   4QJLFOPJTF

    XJUIEPXOTBNQMJOH Y 
  164. Singular value distribution of Hankel matrix 164

  165. In Vivo Motion artifact sudden motion (3 times)

  166. Cardiac Motion artifact

  167. Zipper artifact 167

  168. 2-D herringbone

  169. 2-D herringbone (in-vivo) 169 before after

  170. 2-D herringbone (in-vivo)

  171. Result without k-space weighting 171

  172. EPI Ghost artifact Gx RO Gy PE Gz SS RF

    Ghost artifact image In EPI, Gradient is distorted by eddy currents and this causes phase shift Distorted gradient FT Even and odd echo mismatch causes ghost artifact! Phase shift
  173. Conventional correction •  Navigator : pre-scan or reference scan Navigator-free

    PE RO Make phase dif ference map Navigator-based •  Pulse sequence compensation Calculate difference of phase between 1st -2nd line, 2nd -3rd line only possible to linear phase correction •  Without any modification lower performance compared to the reference-based approaches Gx RO Gy PE Gz SS RF Without PE gradient Xiang QS et al., MRM, 2007 Poser BA et al., MRM, 2013 - Using Parallel Imaging Information Kim YC et al., JMRI, 2007 Zhang et al., MRM, 2004 1) 2) 2) 1) t-1 t … … SENSE recon. SENSE recon. Phase disparity from EPI data itself Calculate - others
  174. EPI model Image intensity Frequency offset Echo time Echo spacing

    (time between each echo) EPI data can be expressed as N : Total # of echoes n : Index of each line x : Read-out y : Phase-encoding Virtual k-space (even signals) Virtual k-space (odd signals) where Different!
  175. Sparsity of difference The ghost generating phase term can be

    changed into a sine term Sparse How can we use this sparsity?
  176. Sparsity of difference (Cont.) Low rank structured matrix completion algorithm

    EPI ghost correction P roblem k-space interpolation Problem using low rank structure
  177. Reconstruction flow •  SE-EPI in-vivo data, 128x128 matrix size, 6/8

    partial Fourier Image K-space Odd echo only Even echo only Reconstruction by odd echo Reconstruction by even echo SSOS result ALOHA
  178. Result : GRE-EPI in-vivo Direct inverse FT Proposed Conventional -with

    reference Conventional -w/o reference Image Re-scaled Image (20%) Phase
  179. Result : fMRI analysis fMRI analysis of GRE-EPI using SPM

    •  Pair hand squeezing stimulation ⇒ Motor cortex activation Proposed Conventional SPMmip [0, -1, 1.25] < < < SPM{T 56 } hand1 SPMresults:.\classical Height threshold T = 5.315495 {p<0.05 (FWE)} Extent threshold k = 0 voxels Design matrix 0.5 1 1.5 2 2.5 10 20 30 40 50 60 contrast 1 Statistics: p-values adjusted for search volume set-level p c cluster-level p FWE-corr q FDR-corr k E p uncorr peak-level p FWE-corr q FDR-corr T (Z ≡ ) p uncorr mm mm mm 0.000 15 0.000 0.000 110 0.000 0.000 0.000 13.22 Inf 0.000 -39 -28 50 0.000 0.000 8.75 6.92 0.000 -45 -13 61 0.000 0.019 6.82 5.80 0.000 -42 -25 65 0.000 0.000 146 0.000 0.000 0.000 13.13 Inf 0.000 45 -19 61 0.000 0.000 9.84 7.46 0.000 45 -22 46 0.000 0.001 7.86 6.43 0.000 51 -10 46 0.000 0.000 36 0.000 0.000 0.000 8.83 6.96 0.000 -51 -37 65 0.000 0.000 11 0.000 0.000 0.019 6.84 5.81 0.000 -6 53 -25 0.000 0.005 6 0.003 0.000 0.019 6.78 5.77 0.000 -39 53 -14 0.000 0.001 10 0.000 0.000 0.025 6.67 5.70 0.000 -51 47 -18 0.000 0.000 13 0.000 0.001 0.088 6.26 5.43 0.000 -3 -4 54 0.000 0.000 11 0.000 0.002 0.097 6.19 5.38 0.000 -39 -52 65 0.021 0.516 5.56 4.94 0.000 -42 -55 58 0.005 0.063 2 0.059 0.003 0.144 6.04 5.28 0.000 -33 59 5 0.005 0.063 2 0.059 0.004 0.154 6.00 5.25 0.000 -3 -94 1 0.005 0.063 2 0.059 0.005 0.188 5.92 5.20 0.000 -54 -31 50 0.005 0.063 2 0.059 0.009 0.285 5.78 5.09 0.000 6 -91 16 Proposed (multi-coil) Conventional Proposed (single-coil)
  180. §  /BOPTDPQZCBTFEPOMPDBMJ[BUJPO §  -PDBMJ[BUJPOQSFDJTJPOJTOPUEJGGSBDUJPOMJNJUFE §  4QBSTFMZBDUJWBUFEQSPCFT MPDBMJ[BUJPOTVQFSS FTPMVUJPOJNBHF §  )PXFWFS

    TQBSTFBDUJWBUJPOTDIFNFIBTUPPTMPXU FNQPSBMSFTPMVUJPOGPSMJWFJNBHJOH §  5FOTPGTFDPOETPSTFWFSBMNJOVUFT §  )JHIEFOTJUZJNBHJOHGPSGBTUMJWFJNBHJOH §  3FRVJSFBSPCVTUMPDBMJ[BUJPOBMHPSJUINBOETZTUFN 6/ Low-density imaging High-density imaging Localization microscopy
  181. Greedy approach Sparsity based approach Min, J.et al, Sci. Rep,

    2014 Zhu, L.et al, Nat Methods, 2012 Holden, S.et al, Nat Methods, 2011 Better Localization Performance 8/ Existing high density algorithm
  182. §  *NBHJOHNPEFM  §  4VDDFTTGVMSFDPWFSJGTPMVUJPOJTTQBSTFCVUEJTDSFUFEPNBJOGPSNVMB Sub-pixel grid c g H

    b Camera image Convolution kernel background Point spread function Sparsity based approach g: measurement, f: sample, c:brightness, b: background, H: point spread function 9/
  183. §  0VSQSFWJPVTBMHPSJUIN'"-$0/ .JOFUBM   §  GBTUMPDBMJ[BUJPOBMHPSJUINCBTFEPODPOUJOVPVTGPSNVMBUJPO §  4NBMMTIJGUFE14'DBOCFBQQSPYJNBUFECZVTJOH5BZMPSTFSJFT 

    §  4QBSTJUZQSPNPUJOHMPDBMJ[BUJPOGPSNVMBPODPOUJOVPVTTQBDF  §  .JOJNJ[BUJPOCZVTJOHGBTUTPMWFSDBMMFE"%.. BMUFSOBUJOHEJSFD UJPONFUIPEPGNVMUJQMJFST  Continuous domain localization 10 /
  184. Fast 2D live cell imaging results §  -JWF&3 FOEPQMBTNJDSFUJDVMVN 1"-.EBUB

    §  5FNQPSBMSFTPMVUJPOTFD GSBNFT  Min et al. Scientific reports (2014) 12 /
  185. §  Hybrid system : Biplane + astigmatism imaging §  PSF

    mutual coherency test •  %.%JDISPJDNJSSPS •  &.&NJTTJPOGJMUFS •  #4#FBN4QMJUUFS •  $-$ZMJOESJDBM-FOT 3D localization and hybrid system ON ON .JOFUBM #JPNFEJDBM0QUJDT&YQSFTT  MFTTBYJBMMZDPIFSFOUCFUUFSBYJBMMPDBMJ[BUJPOQSFDJTJPO  
  186. 3D Reconstruction of fixed cell §  'JYFENJDSPUVCVMF4503.EBUBPGGSBNFT TJOHMFGSBNF XJEFGJFME -FBTUTRVBSFGJUUJOH

    '"-$0/% .JOFUBM #JPNFEJDBM0QUJDT&YQSFTT   
  187. 3D live cell imaging Min et al. Biomedical Optics Express

    (2014) §  -JWF&3 FOEPQMBTNJDSFUJDVMVN 1"-.EBUB §  5IFGBTUFTUMJWFEFNPOTUSBUJPO TIPXJOHUFNQPSBMSFTPMVUJPOTFD 15 /
  188. 188 ALOHA principle ü  14'FTUJNBUJPO ü  %FDPOWPMVUJPO ü  (SJEGSFFMPDBMJ[BUJPO ALOHA

    for localization microscopy (Min et al, TIP, 2018)
  189. 189/ .JOJNVN PSF estimation

  190. 26 / Grid-free localization

  191. 191/ True 3. Grid-free Localization 2. Deconvolution 1. PSF estimation

    Raw data Image Fourier ROI for PSF estimation Algorithm procedure
  192. PSF variation along time

  193. Reconstruction

  194. Localization bias

  195. Summary : Part II •  ALOHA: Measurement domain interpolation using

    the sparsity principle︎ •  ALOHA overcomes the many limitations of CS︎ •  Various biomedical applications of ALOHA ︎ •  Downside: computational complexity︎
  196. References: Part II •  Jong Chul Ye, Jong Min Kim,

    Kyong Hwan Jin and Kiryung Lee, "Compressive sampling using annihilating filter-based low-rank interpolation", IEEE Trans. on Information Theory, vol. 63, no. 2, pp.777-801, Feb. 2017. •  Kyong Hwan Jin, Dongwook Lee, and Jong Chul Ye. "A general framework for compressed sensing and parallel MRI using annihilating filter based low-rank hankel matrix," IEEE Trans. on Computational Imaging, vol 2, no. 4, pp. 480 - 495, Dec. 2016. •  Kyong Hwan Jin, Ji-Yong Um, Dongwook Lee, Juyoung Lee, Sung-Hong Park and Jong Chul Ye, " MRI artifact correction using sparse + low-rank decomposition of annihilating filter-based Hankel matrix", Magnetic Resonance in Medicine (in press), 2016 •  Juyoung Lee, Kyong Hwan Jin, and Jong Chul Ye, "Reference-free single-pass EPI Nyquist ghost correction using annihilating filter-based low rank Hankel matrix (ALOHA)", Magnetic Resonance in Medicine, 10.1002/mrm.26077, Feb. 17, 2016.
  197. References •  Dongwook Lee,, Kyong Hwan Jin, Eung Yeop Kim,

    Sung-Hong Park and Jong Chul Ye, "Acceleration of MR parameter mapping using annihilating filter-based low rank Hankel matrix (ALOHA)", Magnetic Resonance in Medicine, 10.1002/ mrm.26081, Jan. 1, 2016.​ •  Kyong Hwan Jin and Jong Chul Ye, "Annihilating filter based low rank Hankel matrix approach for image inpainting", IEEE Trans. Image Processing, 2015 Nov;24(11):3498-511. •  KH Jin, JC Ye, Sparse+ low rank decomposition of annihilating filter-based Hankel matrix for impulse noise removal, IEEE Trans. on Image Processing, 2018 •  Min, J., Carlini, L., Unser, M., Manley, S., & Ye, J. C. (2015, September). Fast live cell imaging at nanometer scale using annihilating filter-based low-rank Hankel matrix approach. In SPIE Optical Engineering+ Applications (pp. 95970V-95970V). International Society for Optics and Photonics. •  Jin, Kyong Hwan, Yo Seob Han, and Jong Chul Ye. "Compressive dynamic aperture B-mode ultrasound imaging using annihilating filter-based low-rank interpolation." Biomedical Imaging (ISBI), 2016 IEEE 13th International Symposium on. IEEE, 2016.
  198. Jong Chul Ye Bio-Imaging & Signal Processing Lab. Dept. Bio

    & Brain Engineering Dept. Mathematical Sciences KAIST, Korea IEEE EMBS Summer School on Biomedical Imaging Sparse and Deep Learning Approaches for Biomedical Image Reconstruction: Part 3: Deep Learning Updated version of presentation material can be downloaded next week from http://bispl.weebly.com
  199. Compressed Sensing Hankel Structured Matrix Comple7on Deep Learning 199 Roadmap:

    From CS to Deep Learning Coherent Theme of Sparse Recovery
  200. Deep Learning Age •  Deep learning has been successfully used

    for classification, low-level computer vision, etc •  Even outperforms human observers
  201. 201 Reinforcement Learning (AlphaGo)

  202. Generative Adversarial Networks (GoodFellow, 2016) Figure adopted from •  BEGAN:

    Boundary Equilibrium Generative Adversarial Networks, 17’03 •  StackGAN : Text to Photo-realistic Image Synthesis with Stacked Generative Adversarial Networks, ‘16.12
  203. 203

  204. •  Successful demonstraAon of deep learning for various image reconstrucAon

    problems –  Low-dose x-ray CT (Kang et al, Chen et al, Wolterink et al, Ye et al) –  Sparse view CT (Jin et al, Han et al, Adler et al) –  Interior tomography (Han et al) –  Stationary CT for baggage inspection (Han et al) –  CS-MRI (Hammernik et al, Schlemper et al, Yang et al, Lee et al, Zhu et al) –  US imaging (Yoon et al ) –  Diffuse optical tomography (Yoo et al) –  Elastic tomography (Yoo et al) –  Optical diffraction tomography (Kamilov et al) –  etc •  Advantages –  Very fast reconstruction time –  Significantly improved results Deep Learning for Inverse Problems
  205. WavResNet results

  206. WavResNet results

  207. MBIR C D WavResNet results

  208. 1st view 2nd view 3rd view 4th view 5th view

    6th view 7th view 8th view 9th view 9-view Recon
  209. FBP

  210. TV

  211. Deep Learning

  212. WHY DEEP LEARNING WORKS FOR RECON ? DOES IT CREATE

    ANY ARTIFICIAL FEATURES ?
  213. First CNN: LeNet (LeCun, 1998)

  214. Convolutional Layer

  215. Pooling Layer Figure from Leonardo Araujo Aantos

  216. Too Simple to Analyze..? Convolution & pooling à stone age

    tools of signal processing What do they do ?
  217. Dark Age of Applied Mathematics ?

  218. CNN – BIOLOGICAL ORIGIN A LAYMAN’S EXPLANATION

  219. 219 http://klab.smpp.northwestern.edu/wiki/images/4/43/NTM2.pdf Emergence of Hiearchical Features

  220. 220 http://www.vicos.si/File:Lhop-hierarchy-second.jpg •  LeCun et al, Nature, 2015 Hierarchical representation

  221. Visual Information Processing in Brain 221 Kravitz et al, Trends

    in Cognitive Sciences January 2013, Vol. 17, No. 1
  222. Retina, V1 Layer 222 Receptive fields of two ganglion cells

    in retina à convolution Orientation column in V1 http://darioprandi.com/docs/talks/image-reconstruction-recognition/graphics/pinwheels.jpg Figure courtesy by distillery.com
  223. Visual Pathway 223 Poggio et al NATURE|VOL 431 | 14

    OCTOBER 2004
  224. “The Jennifer Anniston Cell” 224 Quiroga et al, Nature, Vol

    435, 24, June 2005
  225. 225

  226. CNN – MATHEMATICAL UNDERSTANDING

  227. Why Deep Learning works for recon ? Existing views 1:

    unfolded iteration •  Most prevailing views •  Direct connecAon to sparse recovery –  Cannot explain the filter channels Jin, arXiv:1611.03679
  228. Why Deep Learning works for recon ? Existing views 2:

    generative model •  Image reconstruc;on as a distribu;on matching –  However, difficult to explain the role of black-box network Bora et al, Compressed Sensing using Generative Models, arXiv:1703.03208
  229. •  What is the role of the nonlinearity such as

    rectified linear unit (ReLU) ? •  Why do we need a pooling and unpooling in some architectures ? •  Why do some networks need fully connected layers whereas the others do not ? •  What is the role of by-pass connection or residual network ? •  What is the role of the filter channels in convolutional layer ? Many Mysteries…
  230. Skipped or Not ? Residual Network Clean image Standard Network

    Zhang, K., et al, IEEE TIP, 2017.
  231. Even Confused about Where to Start Optimization Generalization Architecture

  232. Our Proposal: Deep Learning == Deep Convolutional Framelets •  Ye

    et al, “Deep convolutional framelets: A general deep learning framework for inverse problems”, SIAM Journal Imaging Sciences, 11(2), 991-1048, 2018.
  233. Matrix Representation of CNN Figure courtesy of Shoieb et al,

    2016
  234. Hankel Matrix: Lifting to Higher Dimensional Space

  235. Why we are excited about Hankel matrix ? T -T

    0 n1 -n1 0 * FRI Sampling theory (Ve#erlie et al) and compressed sensing
  236. ︙ ︙ 1 2 3 4 5 6 7 8

    9 -1 0 1 2 3 4 5 6 7 8 9 1 2 3 4 5 6 7 8 9 12 2 3 4 5 6 7 8 9 12 13 3 4 5 6 7 8 9 10 10 10 10 0 11 11 11 1 2 3 4 5 Finite length convolution Matrix Representation * ALOHA : Annihilating filter based LOw rank Hankel matrix Approach * Jin KH et al. IEEE TCI, 2016 * Jin KH et al.,IEEE TIP, 2015 * Ye JC et al. IEEE TIT, 2016 Annihilating filter-based low-rank Hankel matrix
  237. Missing elements can be found by low rank Hankel structured

    matrix compleAon Nuclear norm Projec;on on sampling posi;ons min m kH ( m ) k⇤ subject to P⌦(b) = P⌦( f ) RankH(f) = k * Jin KH et al IEEE TCI, 2016 * Jin KH et al.,IEEE TIP, 2015 * Ye JC et al., IEEE TIT, 2016 m Annihilating filter-based low-rank Hankel matrix
  238. Key Observation Data-Driven Hankel matrix decomposition => Deep Learning • 

    Ye et al, “Deep convolutional framelets: A general deep learning framework for inverse problems”, SIAM Journal Imaging Sciences, 11(2), 991-1048, 2018.
  239. H d(f) = U⌃V T : Non-local basis : Local

    basis Convolution Framelets (Yin et al; 2017) > = I > = I H d(f)
  240. H d(f) H d(f) = ˜ T ˜ T C

    C = T H d(f) C = T (f ~ ) Encoder: ˜ T = I ˜ = PR(V ) H d(f) = U⌃V T Unlifting: f = (˜C) ~ ⌧(˜ ) : Non-local basis : Local basis : Frame condition : rank condition convolution pooling un-pooling convolution : User-defined pooling : Learnable filters H pi (gi) = X k,l [Ci]kl e Bkl i Decoder: Deep Convolutional Framelets (Y, Han, Cha; 2018)
  241. Single Resolution Network Architecture

  242. Multi-Resolution Network Architecture

  243. Network Width versus Depth H d(f) Filter length > Intrinsic

    complexity of signal = Hankel matrix rank Filter length < Intrinsic complexity of signal ? à recursive application of the extended Hankel matrix
  244. Network Width versus Depth

  245. Deep Pyramidal Residual Network Han et al, CVPR, 2017

  246. SO FAR, THE THEORY IS LINEAR.. ROLE OF NONLINEARITY ?

  247. Role of Nonliearity: ReLU à Conic coding D. D. Lee

    and H. S. Seung, Nature, 1999 ReLU: positive framelet coefficients Conic encoding à part by representation similar to visual cortex Han et al, 2018 arXiv:1805.03779
  248. Conic fi [Ci]kl 0 H pi (gi) = X k,l

    [Ci]kl e Bkl i H pi (fi) ' Lifting Geometry of Neural Network gi Un-lifting Ci(fi) Ci(fi) 0 i ⇣i ⇣ e i ⌘
  249. fi [Ci]kl 0 H pi (gi) = X k,l [Ci]kl

    e Bkl i H pi (fi) ' Lifting Geometry of Residual Neural Network Ci(fi) Ci(fi) 0 i ⇣i ⇣ e i ⌘ gi Un-lifting
  250. https://www.youtube.com/watch?v=DdC0QN6f3G4 Relativity: Lifting to the Space-time à linear trajectory High

    Dimensional Lifting in Relativity Falling apple Apple at rest Universal law of gravity: 3-D Space à curved trajectory Falling apple Apple at rest
  251. Deep CNN Lifting Un-lifting Conic Lifting Un-lifting Lifting Un-lifting

  252. Successive Low-Rank Hankel Matrix Approximation Truncated channel à Low rank

    Hankel matrix approximation
  253. Compressed Sensing Hankel Structured Matrix Comple7on Deep Learning 253 From

    CS to Deep Learning: Coherent Theme
  254. APPLICATION-DRIVEN EVIDENCES

  255. U-Net Architecture (Ronnenberger, 2015)

  256. Problem of U-net Pooling does NOT satisfy the frame condition

    JC Ye et al, SIAM Journal Imaging Sciences, 2018 Y. Han et al, TMI, 2018. ext > ext = I + > 6= I
  257. Improving U-net using Deep Conv Framelets •  Dual Frame U-net

    •  Tight Frame U-net JC Ye et al, SIAM Journal Imaging Sciences, 2018 Y. Han and J. C. Ye, TMI, 2018
  258. U-Net versus Dual Frame U-Net

  259. Tight-Frame U-Net JC Ye et al, SIAM Journal Imaging Sciences,

    2018
  260. Denoising: U-Net vs. Tight-Frame U-Net

  261. Inpainting: U-Net vs. Tight-Frame U-Net

  262. k-Space Deep Learning for Accelerated MRI Han et al, arXiv:1805.03779

    ALOHA : k-space interpolation à k-space interpolation using deep learning ? Yes
  263. k-Space Deep Learning for Accelerated MRI Han et al, arXiv:1805.03779

    ~3dB gain
  264. k-Space Deep Learning for Parallel MRI Han et al, arXiv:1805.03779

  265. k-Space Deep Learning for Parallel MRI Cha et al, arXiv:1806.00806

    ~13dB gain
  266. k-Space Deep Learning for Parallel MRI Cha et al, arXiv:1806.00806

  267. DEEP NETWORKS FOR CT

  268. Low-Dose CT •  To reduce the radiation exposure, sparse-view CT,

    low-dose CT and interior tomography. Sparse-view CT (Down-sampled View) Low-dose CT (Reduced X-ray dose) Interior Tomography (Truncated FOV)
  269. 269 Wavelet transform level 2 level 1 level 3 level

    4 Wavelet recomposition + Residual learning : Low-resolution image bypass High SNR band CNN (Kang, et al, Medical Physics 44(10))
  270. Routine dose Quarter dose (Kang, et al, Medical Physics 44(10)

    2017)
  271. Routine dose AAPM-Net results (Kang, et al, Medical Physics 44(10)

    2017)
  272. WavResNet results (Kang et al, TMI, 2018)

  273. WavResNet results (Kang et al, TMI, 2018)

  274. MBIR Our latest Result C D WavResNet results

  275. Full dose Quarter dose

  276. Full dose Quarter dose

  277. Sparse-View CT •  To reduce the radiation exposure, sparse-view CT,

    low-dose CT and interior tomography. Sparse-view CT (Down-sampled View) Low-dose CT (Reduced X-ray dose) Interior Tomography (Truncated FOV)
  278. Multi-resolution Analysis & Receptive Fields

  279. Tight-Frame U-Net JC Ye et al, SIAM Journal Imaging Sciences,

    2018 Han et al, TMI, 2018
  280. 90 view recon U-Net vs. Tight-Frame U-Net •  JC Ye

    et al, SIAM Journal Imaging Sciences, 2018 •  Y. Han and J. C. Ye, TMI, 2018
  281. None
  282. None
  283. •  Figures from internet 9 View CT for Baggage Screening

  284. 9 View CT for Baggage Screening

  285. 9 View Dual Energy CT for Baggage Screening Han et

    al, arXiv preprint arXiv:1712.10248, (2017). CT meeting 2018
  286. 1st view 2nd view 3rd view 4th view 5th view

    6th view 7th view 8th view 9th view
  287. FBP

  288. TV

  289. Ours

  290. ROI Reconstruction •  To reduce the radiation exposure, sparse-view CT,

    low-dose CT and interior tomography. Sparse-view CT (Down-sampled View) Low-dose CT (Reduced X-ray dose) Interior Tomography (Truncated FOV)
  291. Deep Learning Interior Tomography Han et al, arXiv preprint arXiv:1712.10248,

    (2017): CT meeting 2018.
  292. Ground Truth FBP

  293. TV Chord Line Ours 8~10 dB gain

  294. Ground Truth Chord Line TV Ours

  295. Summary: Part III •  Deep learning for inverse problems︎ • 

    Significant performance gain︎ •  Deep convolutional framelets: new mathematical tools for understanding deep neural network for inverse problems︎ •  Depth vs. width︎ •  Role of ReLU︎ •  Residual net︎ •  Numerical results support the theory︎
  296. Still Unresolved Problems.. •  Cascaded geometry of deep neural network

    •  Generalization capability •  Optimization landscape •  Training procedure •  Extension to classification problems
  297. LAB CONTENTS

  298. Dataset of Natural Images •  MNIST (http://yann.lecun.com/exdb/mnist/) –  Handwritten digits.

    •  SVHN (http://ufldl.stanford.edu/housenumbers/) –  House numbers from Google Street View. •  ImageNet (http://www.image-net.org/) –  The de-facto image dataset. •  LSUN (http://lsun.cs.princeton.edu/2016/) –  Large-scale scene understanding challenge. •  Pascal VOC ( http://host.robots.ox.ac.uk/pascal/VOC/) –  Standardized image dataset. •  MS COCO (http://cocodataset.org/#home) –  Common Objects in Context. •  CIFAR-10 / -100 (https://www.cs.utoronto.ca/~kriz/cifar.html) –  Tiny images data set. •  BSDS300 / 500 ( https://www2.eecs.berkeley.edu/Research/ Projects/CS/vision/grouping/resources.html) –  Contour detection and image Segmentation resources. https://en.wikipedia.org/wiki/List_of_datasets_for_machine_learning_research https://www.kaggle.com/datasets
  299. Dataset for Medical Images •  HCP (https://www.humanconnectome.org/) –  Behavioral and

    3T / 7T MR imaging dataset. •  MRI Data (http://mridata.org/) –  Raw k-space dataset acquired on a GE clinical 3T scanner. •  LUNA (https://luna16.grand-challenge.org/data/) –  Lung Nodule analysis dataset acquired on a CT scanner. •  Data Science Bowl (https://www.kaggle.com/c/data-science-bowl-2017) –  A thousand low-dose CT images. •  NIH Chest X-rays (https://nihcc.app.box.com/v/ChestXray-NIHCC) –  X-ray images with disease labels. http://www.cancerimagingarchive.net/ https://www.kaggle.com/datasets •  TCIA Collections (http://www.cancerimagingarchive.net/) •  De-identifies and hosts a large archive of medical images of cancer accessible for public download. •  The data are organized as “Collections”, typically patients related by a common disease, image modality (MRI, CT, etc).
  300. Libraries for Deep learning •  TensorFlow (https://www.tensorflow.org/) –  Python • 

    Theano ( http://deeplearning.net/software/theano/) –  Python •  Keras (https://keras.io/) –  Python •  Caffe (http://caffe.berkeleyvision.org/) –  Python •  Torch (or PyTorch) (http://torch.ch/) –  C / C++ (or Python) •  Deeplearning4J ( https://deeplearning4j.org/) –  Java •  Microsoft Cognitive Toolkit (CNTK) ( https://www.microsoft.com/en-us/cognitive- toolkit/) –  Python / C / C++ •  MatConvNet ( http://www.vlfeat.org/matconvnet/) –  Matlab
  301. Implementation on MatConvNet Download the MatConvNet Provide the pre-trained models

  302. Step 1: Compile the toolbox 1.  Unzip the MatConvNet toolbox.

    2.  Open ‘vl_compilenn.m’ in Matlab.
  303. Step 1: Compile the toolbox (cont.) 3.  Check the options

    such as enableGpu and enableCudnn. 4.  Run the ‘vl_compilenn.m’. * To use GPU processing (false true), you must have CUDA installed. ( https://developer.nvidia.com/cuda-90-download-archive ) ** To use cuDNN library (false true), you must have cuDNN installed. ( https://developer.nvidia.com/cudnn )
  304. Step 2: Prepare dataset As a classification example, MNIST consists

    of as follows, Images % struct-type Data % handwritten digit image labels % [1, …, 10] set % [1, 2, 3], 1, 2, and 3 indicate train, valid, and test set, respectively. Data Labels 6
  305. Step 2: Prepare dataset (cont.) •  As a segmentation example,

    U-Net dataset consists of as follows, –  Images % struct-type Ø  data % Cell image Ø  labels % [1, 2], Mask image. 1 and 2 indicate back- and for-ground, respectively. Ø  set % [1, 2, 3] Data Labels
  306. Step 3: Implementation of the network architecture •  Developers only

    need to program the network architecture code because MatConvNet supports the network training framework. Support famous network architectures, such as alexnet, vggnet, resnet, inceptionent, and so on.
  307. Step 3: Implementation of the architecture (cont.) –  The implementation

    details of U-Net U-Net can be implemented, recursively. Stage 0 Stage 1 Stage 2 Stage 3 Stage 4
  308. Step 3: Implementation of the architecture (cont.) 1.  Create objects

    of network and layers. Encoder Part Skip + Concat Part Decoder Part •  The structure of Stage 0 Network Part
  309. Step 3: Implementation of the architecture (cont.) 2.  Connect each

    layers. •  The structure of Stage 0 Layer Name ( string-type ) Layer object ( object ) Input Name ( string-type ) Output Name ( string-type ) Parameters Name ( string-type ) All objects and names must be unique.
  310. Step 3: Implementation of the architecture (cont.) 3.  Implement recurrently

    the each stages and add a loss function. Previous parts (3.1 and 3.2) become functional as ‘add_block_unet’.
  311. Step 4: Network hyper-parameter set up •  MatConvNet supports the

    default hyper-parameters as follows, Refer the cnn_train.m ( or cnn_train_dag.m ) The supported hyper-parameters 1.  The size of mini-batch 2.  The number of epochs 3.  Learning rate 4.  Weight decay factor 5.  Solvers such as SGD, AdaDelta, AdaGrad, Adam, and RMS The kind of Optimization Solvers
  312. Step 5: Run the training script 1. Training script 2.

    Training loss 3. Training loss graph •  Blue : train •  Orange : valid
  313. Cell segmentation

  314. train_cnn_cell_seg.m val_cnn_cell_seg.m

  315. F5 ژח Run // Training code Training code

  316. Ѿҗ : figure 1

  317. Ѿҗ : figure 10

  318. F5 ژח Run ߡౡ ௿ܼ // Validation code Validation code

    // ޷ܻ ੷੢೧֬਷ ֎௼ਕ௼ܳ ࠛ۞৬ प೯
  319. F5 ژח Run ߡౡ ௿ܼ Validation code // ೟णػ ݃૑݄

    ֎௼ਕ௼ܳ ࠛ۞৬ प೯
  320. Ѿҗ : figure 20

  321. References: Part III •  Ge Wang, Jong Chu Ye, Klaus

    Mueller, Jeffrey A Fessler, "Image Reconstruction Is a New Frontier of Machine Learning", IEEE Trans. on Medical Imaging, Vol. 37 no. 6, pp. 1289 - 1296, June 2018. •  Yoseob Han and Jong Chul Ye,"Framing U-Net via Deep Convolutional Framelets: Application to Sparse-view CT", Special Issue on Machine Learning for Image Reconstruction, IEEE Trans. on Medical Imaging, vol. 37, no. 6, pp. 1418-1429, June, 2018. •  Eunhee Kang, Won Chang, Jaejun Yoo, and Jong Chul Ye,"Deep Convolutional Framelet Denosing for Low-Dose CT via Wavelet Residual Network", Special Issue on Machine Learning for Image Reconstruction, IEEE Trans. on Medical Imaging, vol. 37, no.6, pp. 1358-1369, 2018. •  Dongwook Lee, Jaejun Yoo, Sungho Tak and Jong Chul Ye, "Deep Residual Learning for Accelerated MRI using Magnitude and Phase Networks", IEEE Trans on Biomedical Engineering (in press), Invited paper for Special Section on Deep Learning, 2018. •  Jong Chul Ye, Yoseob Han and Eunju Cha, "Deep convolutional framelets: a general deep learning framework for inverse problems", SIAM Journal on Imaging Sciences (in press), 2018. •  Yoseob Han, Jaejun Yoo, Hak Hee Kim, Hee Jung Shin, Kyunghyun Sung, and Jong Chul Ye, "Deep Learning with Domain Adaptation for Accelerated Projection-Reconstruction MR", Magnetic Resonance in Medicine (accepted), 2018. •  Eunhee Kang, Junhong Min and Jong Chul Ye, " A Deep Convolutional Neural Network using Directional Wavelets for Low-dose X-ray CT Reconstruction", Medical Physics 44.10 (2017).
  322. Acknowledgements CT Team •  Yoseob Han •  Eunhee Kang • 

    Jawook Goo US Team •  Shujaat Khan •  Jaeyong Hur MR Team •  Dongwook Lee •  Juyoung Lee •  Eunju Cha •  Byung-hoon Kim Image Analysis Team •  Boa Kim •  Junyoung Kim Optics Team •  Sungjun Lim •  Junyoung Kim •  Jungsol Kim •  Taesung Kwon