Lecture 5: Transforms, Fourier and
Wavelets
Outline
Talk involves matrices and vector spaces
– You will not be tested on it
What are Transforms
= change of basis
Linear or non-linear (will focus on linear)
Fourier Transforms
Wavelet Transforms
Why??
– Because a transform or a change in basis may
allow you to see things differently, see things that
couldn’t be seen before, to get a different
“perspective”
IDEA Lab, Radiology, Cornell 2
What are transforms?
Vectorize a signal (ECG, MR image, …) into vector x
A linear transform on this vector is defined as a matrix operation
y = Tx
– Linearity: T(x1 + x2) = T x1 + T x2
Matrix examples
T is generally a square, full-rank matrix
If T is a “wide” matrix, then the transform does not have a
unique inverse
– Also known as overcomplete transform
• T is orthogonal if Tt T = diagonal matrix
• T is orthonormal if Tt T = Identity matrix
• Orthonormal transforms retain signal energy
||Tx|| = ||x||
IDEA Lab, Radiology, Cornell 3
Transforms
Examples:
– Fourier transform is an orthonormal transform
– Wavelet transform is generally overcomplete, but
there also exist orthonormal wavelet transforms
A good property of a transform is invertibility
– Both Fourier and wavelet transforms are invertible
Many other image-based processes are not
invertible
– E.g. Distance transform, JPEG compression, edge
detection, blurring
IDEA Lab, Radiology, Cornell 4
Transforms
A transform (with full rank T) is a change of basis
Definition: A basis on a vector space is a set of
linearly independent vectors that are able to
express any other vector of the space as a linear
combination of them.
v
c2 e2 v = c1 e1 + c2 e2
c1 e1
5
IDEA Lab, Radiology, Cornell
Basis and basis vectors
A set of such vectors can be expressed as a matrix
T = [t1 | t2 | … tn], where each tk is a column vector
Basis Basis
matrix vectors
Any
Any subset
subset of
of the
the vectors
vectors T
TKK =
= {t
{tkk}
}k=1:K then defines a
k=1:K then defines a
subspace
subspace of Rnn
of R
–– Technically
Technically wewe say
say the
the subspace
subspace is is spanned
spanned byby the
the
vectors
vectors {t
{tkk}
} or
or equivalently
equivalently byby the
the matrix
matrix T
TKK
Many
Many basis
basis matrices
matrices can
can span
span the
the same
same space
space
IDEA Lab, Radiology, Cornell 6
Basis
Example in R2:
– a pair of vectors rotated by 45 degree:
1 / 2 −1/ 2 1/ 2 1/ 2
R = R −1 =
1/ 2
1 / 2 1 / 2
−1/ 2 1/ 2
1/ 2
– If you want to change any vector in standard space to this new
basis, just left-multiply it by R-1
– Note: R is orthonormal: R-1 = Rt
IDEA Lab, Radiology, Cornell 7
Transform = change of basis
T is a basis if its columns are linearly independent
– i.e. it is a full rank matrix
Orthonormality simply makes everything easier
T-1 = Tt or for complex vectors, T-1 = TH
So computing the inverse transform is very easy:
If y = Tx, then x = TH y
This is why search for (useful) orthonormal
transforms is such a huge deal
IDEA Lab, Radiology, Cornell 8
Other examples in medical imaging
Radon transform widely used to turn raw CT
data into CT images
– X-ray absorption is a line integral
Funk-Radon is an extension of it, and is used to
reconstruct orientation distribution function
(ODF) from diffusion MRI data
Another transform (spherical harmonic
transform) is used to clean up ODF
IDEA Lab, Radiology, Cornell 9
High Angular Resolution Diffusion Imaging
IDEA Lab, Radiology, Cornell 10
MR Diffusion Imaging
Diffusion MRI measures the directionally varying diffusion
properties of water in tissue
D-MRI involves taking several directional diffusion
imaging measurements
Then we fit a 3D shape to these measurements
IDEA Lab, Radiology, Cornell 11
•Data Acquisition Strategy
•RF
•Gs
•qz
•qy
•Gx
•qx
•Gy
IDEA Lab, Radiology, Cornell
Reconstruction Problem
Basic Approach Construct a function on the unit sphere
characterizing the angular structure of diffusion in each voxel.
•S(φ,θ •F(φ,θ
) )
Recon using spherical harmonic basis Let f and s be vectors
representing functions S(.) and F(.). Then
f = [SH transform] [F-R transform] s
Represents ODF as linear
mix of spherical harmonics Transforms raw MR data to
function on unit sphere
IDEA Lab, Radiology, Cornell
•High Angular Resolution Diffusion Imaging:
Spherical Harmonic Transform
Point Spread Functions Example ODFs
= Basis functions or vectors
A linear combination of these
(and their rotated versions)
t1 t2 t3 t4
can construct an arbitrary
ODF on the unit sphere
IDEA Lab, Radiology, Cornell
• Middle cerebellar peduncle (MCP)
• Superior cerebellar peduncle (SCP)
• Pyramidal tract (PT)
• Trans pontocerebellar fibers (TPF)
Hess CP, Mukherjee P, Han ET, Xu D, Vigneron DB
IDEA Lab, Radiology,
Magn ResonCornell
Med 2006; 56:104-117 15
Fourier Transforms - Audio example
Audio signals like music have various frequencies
– You can change the contribution of various
frequency bands by using a band equalizer
Each frequency band is represented by a pure
tone at frequency k Hz
Lets say its captured in a vector ek
The contribution of this band = the dot product
sk = <ek, x> = ekH x
IDEA Lab, Radiology, Cornell 16
Fourier Transforms (FT)
Now do this for each k, and collect the sk‘s in vector s, we have:
s = FH x
This is called the Fourier Transform
F = (e1 | e 2 | ... e n )
i 2π
Basis vectors exp( − k .1)
n Complex exponentials (like sinusoids)
where e k = Μ
exp( − i 2π k .n )
n
Class task: show that F is an orthonormal basis
Therefore, FT is easily invertible:
x=FX
FT = change of basis from standard space to Fourier space
IDEA Lab, Radiology, Cornell 17
What is Fourier space and why
does anyone need it?
Fourier basis is a collection of harmonics
– Note that complex exponentials are simply sines and
cosines
Therefore the FT simply decomposes a signal into its
harmonic components
FT gives direct information about the sharpness and
oscillations present in the data
An “alternate view” of the data
IDEA Lab, Radiology, Cornell 18
FT Demo
See how a unit pulse signal is constructed from
Fourier basis vectors
Fourier coefficients
Basis vectors
Reconstructed signal
IDEA Lab, Radiology, Cornell 19
Fast Fourier Transform
The computation of FT requires n2 mult operations
Can get pretty expensive for large n
An efficient algorithm exists which can do the job in
O(n log(n)) operations
This is extremely fast vs original FT for large n
Most programming languages have a FFT library
– In C++ use FFTW, in MATLAB, built-in function fft.m
IDEA Lab, Radiology, Cornell 20
FT in images
FT is defined on 1D, 2D or nD data.
In 2D for instance you do FT along image rows,
then do FT along columns
Again, the FT coefficients are dot products of the
image with complex exponential basis vectors
Basis vectors represent frequency (spatial
frequency, or how “sharp things are)
The FT coefficients represent the contribution of
each spatial frequency
Fourier space in 2D is sometimes called “k-space”
IDEA Lab, Radiology, Cornell 21
k-space and image-space
k-space &
image-
space are
related by
the 2D FT
IDEA Lab, Radiology, Cornell 22
Truncation
Just as the number of frequency bands
determines the highest pitch in an audio signal,
the number of k-space points determines the
sharpest features of an image
Truncation = sampling central part of k-space
IDEA Lab, Radiology, Cornell 23
Ringing Example
IDEA Lab, Radiology, Cornell 24
Further Reading
IDEA Lab, Radiology, Cornell 25
Wavelet Transform
FTs are great, but they capture global features
– Harmonic components of the entire signal
– They are obtained by dot-producting the WHOLE signal
Problem1: local features can get lost
Problem2: if signal is not stationary (features change
with time or in space) then this is not captured by FT
Therefore need a transform that provides frequency
information LOCALLY
IDEA Lab, Radiology, Cornell 26
Wavelet Transforms - motivation
For example consider the following signal
x(t)=cos(2*pi*10*t)+cos(2*pi*25*t)+cos(2*pi*50*t)+cos(2*pi*100*t)
Has frequencies 10, 25, 50, and 100 Hz at any
given time instant
stationary signal, FT
can provide full info
Non-stationary, frequency
content changes with time
FT CANNOT provide full info
IDEA Lab, Radiology, Cornell 27
Time-frequency, time-scale analysis
What we need is a time-frequency analysis
Do FT in a local time window
generalization of
local frequency
Note different tiling
frequency scale
of t-s space
Insight: time window
depends on scale
time time
IDEA Lab, Radiology, Cornell 28
Basis functions in WT
Haar
Mexican Hat
Daubechies
(orthonormal)
IDEA Lab, Radiology, Cornell 29
WT in images
Images are piecewise smooth or piecewise constant
Stationarity is even rarer than in 1D signals
FT even less useful (nnd WT more attractive)
2D wavelet transforms are simple extensions of 1D WT,
generally performing 1D WT along rows, then columns etc
Sometimes we use 2D wavelets directly, e.g. orthonormal
Daubechies 2D wavelet
IDEA Lab, Radiology, Cornell 30
WT on images
2D generalization of scale-
Scale 0
time decomposition
H
scale Scale 1
time Scale 2
V H-V
Successive application of dot product with wavelet of increasing width.
Forms a natural pyramid structure. At each scale:
H = dot product of image rows with wavelet
V = dot product of image rows with wavelet
H-V = dot product of image rows then columns with wavelet
IDEA Lab, Radiology, Cornell 31
Wavelet Applications
Many, many applications!
Audio, image and video compression
New JPEG standard includes wavelet compression
FBI’s fingerprints database saved as wavelet-
compressed
Signal denoising, interpolation, image zooming,
texture analysis, time-scale feature extraction
In our context, WT will be used primarily as a
feature extraction tool
Remember, WT is just a change of basis, in order
to extract useful information which might
otherwise not be easily seen
IDEA Lab, Radiology, Cornell 32
WT in MATLAB
MATLAB has an extensive wavelet toolbox
Type help wavelet in MATLAB command window
Look at their wavelet demo
Play with Haar, Mexican hat and Daubechies wavelets
IDEA Lab, Radiology, Cornell 33
Project Ideas
Idea 1: use WT to extract features from ECG data
– use these features for classification
Idea 2: use 2D WT to extract spatio-temporal
features from 3D+time MRI data
– to detect tumors / classify benign vs malignant tumors
Idea 3: use 2D WT to denoise a given image
IDEA Lab, Radiology, Cornell 34
Idea 3: Voxel labeling from
contrast-enhanced MRI
Can segment according to time profile of 3D+time
contrast enhanced MR data of liver / mammography
Temporal models used to
extract features
Typical plot of time-resolved
MR signal of various tissue Instead of such a simple temporal model,
classes wavelet decomposition could provide
spatio-temporal features that you can
use for clustering
IDEA Lab, Radiology, Cornell
Liver tumour quantification from
DCE-MRI
IDEA Lab, Radiology, Cornell
Further Reading on Wavelets
IDEA Lab, Radiology, Cornell 37
Lecture 5: Transforms, Fourier and
Wavelets
