# Matlab - PowerPoint

Document Sample

Matlab and Medical Image
Analysis basics / Human Image
Perception

Kostas Marias
Today‘s goals

• Learn enough matlab to get started.
• Review some basics of variables and
algebra in general
• Run some code together
Introduction
What is MATLAB?
• MATLAB is a tool for doing numerical
computations with matrices and vectors. It is
very powerful and easy to use. In fact, it
integrates computation, visualization and
programming all together in an easy-to-use
environment and can be used on almost all the
platforms: windows, Unix, and Apple Macintosh,
etc.

• MATLAB stands for ―Matrix Laboratory‖.
What is Matlab?
• A software environment for interactive numerical computations

• Examples:
– Matrix computations and linear algebra
– Solving nonlinear equations
– Numerical solution of differential equations
– Mathematical optimization
– Statistics and data analysis
– Signal processing
– Modelling of dynamical systems
– Solving partial differential equations
– Simulation of engineering systems
Matlab Background
•Matlab = Matrix Laboratory
•Originally a user interface for numerical linear algebra routines
(Lapak/Linpak)
•Commercialized 1984 by The Mathworks
•Since then heavily extended (defacto-standard)

•Alternatives                             Complements
•Matrix-X                                 Maple            (symbolic)
Octave           (free; GNU)              Mathematica      (symbolic)
Lyme             (free; Palm)
MATLAB
Environment

To start MATLAB:
START  PROGRAMS
 MATLAB 7.0 
MATLAB 7.0
Display Windows
Help!
• help
• help command
Eg., help plot
• Help on toolbar
• demo
• http://math.ucsd.edu/~driver/21d-s99/matlab-primer.html
Matlab Workspace
• Variables are stored here in the form of numerical
matrices which can be of different types: int8,
uint8, int16, uint16, double.
– All numerical operations need to be performed
on double precision, takes 8 bytes per value.
– Other types are efficient for storage (<= 2 bytes)
– Who, whos – current variables in workspace
– Clear all – clear all variables
Matrices in Matlab
• How to build a matrix?
– >>A=[1 2 3; 4 5 6; 7 8 9]; 3x3 matrix
• Special matrices:
– zeros(r,c), ones(r,c), eye(r,c), rand(r,c)
…
• Acess matrix elements:
– A(2,3) - result = 6
– using one number to index: A(8) – result = 6 –
matrix is written out as a single column.
– A(:,2) – obtain second column of the matrix
Matrices in Matlab (2)
• Basic Operations defined on matrices:
– +,-,*,/,^,’,sqrt,sin,cos, etc.
– Element-wise, precede with dot: .*, ./, .^
etc.
• size(A) – size vector
• sum(A) – columns sum vector
• sum(sum(A)) – sum of all elements
• Type help elmat for more instructions
Programming in Matlab
• Script Files
– List of commands to be executed sequentially. Useful
when same sequence is executed many times. They
should be saved with extension script.m and should
be placed on work path to be called.
• >>script   runs the sequence
• Functions
– Similar to script, but takes arguments
– A funcname.m file starts with:
function [output]=funcname(input)
– Any function can be listed with: >>type funcname
• See Matlab editor (File>New>m-file)
Interactive Calculations

• Matlab is interactive, no need to declare variables
• >> 2+3*4/2
• >> a=5e-3; b=1; a+b

•   Most elementary functions and constants are already defined
•   >> cos(pi)
•   >> abs(1+i)
•   >> sin(pi)

• Last call gives answer 1.2246e-016 !?
What kind of graphics is possible in Matlab?

Polar plot:
Line plot:
t=0:.01:2*pi;
polar(t,abs(sin(2*t).*cos(2*t)));         x=0:0.05:5;,y=sin(x.^2);,plot(x,y);

Stem plot:
x = 0:0.1:4;, y = sin(x.^2).*exp(-x);
stem(x,y)
What kind of graphics is possible in Matlab?

Mesh plot:                 Surface plot:
z=peaks(25);, mesh(z);     z=peaks(25);, surf(z);, colormap(jet);

Contour plot:
Quiver plot:
z=peaks(25);,contour(z,16);
Vectors and Matrices

• Vectors (arrays) are defined as
• >> v = [1, 2, 4, 5]
• >> w = [1; 2; 4; 5]

• Matrices (2D arrays) defined similarly
• >> A = [1,2,3;4,-5,6;5,-6,7]
Matrix Operators
• All common operators are overloaded
• >> v + 2

•   Common operators are available
•   >> B = A’
•   >> A*B
•   >> A+B

• Note:
• Matlab is case-sensitive
A and a are two different variables
• Transponate conjugates complex entries; avoided
by
•     >> B=A.’
Indexing Matrices
•Indexing using parentheses
•>> A(2,3)

•Index submatrices using vectors
of row and column indices
•>> A([2 3],[1 2])

•Ordering of indices is important!
•>> B=A([3 2],[2 1])
•>> B=[A(3,2),A(3,1);A(2,2);A(2,1)]
Indexing Matrices
•Index complete row or column using
the colon operator
•>> A(1,:)

•Can also add limit index range
•>> A(1:2,:)
•>> A([1 2],:)

•General notation for colon operator
•>> v=1:5
•>> w=1:2:5
Matrix Functions
• Many elementary matrices predefined
• >> help elmat;
• >> I=eye(3)

• Elementary functions are often overloaded
• >> help elmat
• >> sin(A)

•   Specialized matrix functions and operators
•   >> As=sqrtm(A)
•   >> As^2
•   >> A.*A

• Note: in general, ‖.<operator>‖ is elementwise
operation
Numerical Linear Algebra
• Basic numerical linear algebra
• >> z=[1;2;3]; x=inv(A)*z
• >> x=A\z

•   Many standard functions predefined
•   >> det(A)
•   >> rank(A)
•   >> eig(A)

• The number of input/output arguments can often be
varied
• >> [V,D]=eig(A)
Graphics
•   Visualization of vector data is available
•   >> x=-pi:0.1:pi; y=sin(x);
•   >> plot(x,y)
•   >> plot(x,y,’s-’)
•   >> xlabel(’x’); ylabel(’y=sin(x)’);

• Can change plot properties in Figure menu, or
via ‖handle‖
• >> h=plot(x,y); set(h, ’LineWidth’,
4);
Graphics
•Three-dimensional graphics
•>> A = zeros(32);
•>> A(14:16,14:16) = ones(3);
•>> F=abs(fft2(A));
•>> mesh(F)
•>> rotate3d on

•Several other plot functions available
•>> surfl(F)

•Can change lightning and material
properties
•>> material metal
Graphics
•   Bitmap images can also be visualized
•   >> image(X); colormap(map)
•   >> axis image off
Basic concepts from
algebra
`n
 v

`n
 k
`n
 3

`n
 1          `n
 2
Inner (dot) Product
v     
w        v.w  ( x1 , x2 ).( y1 , y2 )  x1 y1  x2 . y2

The inner product is a SCALAR!

v.w  ( x1 , x2 ).( y1 , y2 ) || v ||  || w || cos

v.w  0  v  w
Matrices
 a11 a12    a1m 
a           a2 m 
Sum:
 21 a22                     Cnm  Anm  Bnm
Anm    a31 a32     a3m 
                  
                         cij  aij  bij
an1 an 2
            anm 
A and B must have the same
dimensions
Matrices
Product:                   A and B must have
Cn p  Anm Bm p   compatible dimensions
m
cij   aik bkj
k 1         Ann Bnn  Bnn Ann

1 0  0
Identity Matrix:
         
 0 1  0  IA  AI  A
I 
   

0 0  1 
         
Matrices
Transpose:
Cmn  A nm
T
( A  B)  A  B
T    T     T

cij  a ji            ( AB)T  BT AT

If    AT  A      A is symmetric
Matrices
Determinant:         A must be square

 a11 a12  a11 a12
det           a         a11a22  a21a12
a21 a22      21 a22

 a11 a12   a13 
det a21 a22          a a22   a23           a21 a23           a21 a22
           a23  11               a12              a13
a32   a33           a31 a33           a31 a32
a31 a32
               
a33 
Matrices
Inverse:                                         A must be square
1              1
Ann A        nn   A        nn   Ann  I
1
 a11 a12           1                 a22  a12 
a          a a a a                 a        
 21 a22      11 22   21 12           21 a11 
Matrix

MATLAB works with essentially only one
kind of object – a rectangular numerical
matrix with possible complex entries.
Entering a Matrix
Matrices can be
 Entered manually;

 Generated by built-in
functions;

An Example
• A = [1, 2, 3; 7, 8, 9]

§        Use ‗ ; ‘ to indicate the end of each
row

§     Use comma or space to separate
elements of a row
Matrix operations:
•   - subtraction
•   * multiplication
•   ^ power
•   ‗ transpose
•   \ left division, / division
x = A \ b is the solution of A * x = b
x = b / A is the solution of x * A = b

• To make the ‗*‘ , ‗^‘, ‗\‘ and ‗/‘ entry-
wise, we precede the operators by ‗.‘
Example:

A = [1, 2, 3; 4, 5, 6]
B = [0.5; 1;3]

   C = A * B v.s.
   C = A .* B
Array Operations
•   Scalar-Array Mathematics
For addition, subtraction, multiplication, and division of an array by a
scalar simply apply the operations to all elements of the array.
•   Example:
>> f = [ 1 2; 3 4]
f=
1 2
3 4
>> g = 2*f – 1
g=
Each element in the array f is
1     3              multiplied by 2, then subtracted
5     7              by 1.
Array Operations (con‘t…)
• Element-by-Element Array-Array Mathematics.

Operation                 Algebraic Form    MATLAB

Subtraction                   a–b            a–b
Multiplication                axb            a .* b
Division                      ab            a ./ b
Exponentiation                  ab           a .^ b

• Example:
Each element in x is multiplied by
>> x = [ 1 2 3 ];
>> y = [ 4 5 6 ];
>> z = x .* y        the corresponding element in y.
z=
4 10 18
Algebraic operations in Matlab:

Array products: Recall that addition and subtraction of matrices involved addition or
subtraction of the individual elements of the matrices. Sometimes it is desired to simply
multiply or divide each element of an matrix by the corresponding element of another
matrix 'array operations”.

Array or element-by-element operations are executed when the operator is preceded
by a '.' (period):

a .* b   multiplies each element of a by the respective element of b
a ./ b   divides each element of a by the respective element of b
a .\ b   divides each element of b by the respective element of a
a .^ b   raise each element of a by the respective b element
Algebraic operations in Matlab:

For example, if matrices G and H are G = [ 1 3 5; 2 4 6]; H = [-4 0 3; 1 9 8];

D=G .* H = [ -4 0 15
2 36 48 ]
Subscripts:

• Subscripts: the element in row i and
column j of A is denoted by A(i, j).

• Example: A = magic(3);
A(1,1) + A(1,2) + A(1,3)
The Colon Operator ‗:‘
The colon ‗:‘ is one of MATLAB ‘s most
important operators. It has many formats:

• 3 : -2 : -11 is a row vector containing integers
from 3 to –11 with a increment of –2.
• Subscript expressions involving colons refer to
portions of a matrix:
A(1:3 , 2) is the first to the third elements of the
second column of A.

What does A([1:3], [2:4]) mean?
Working with Matrices:
MATLAB provides four functions that generate
basic matrices:

• Zeros: all zeros.     A = zeros(1,3)
• Ones: all ones. A = ones(2,4)
• Rand: uniformly distributed random elements.
A = rand(3,5)
• Randn: normally distributed random elements.
A = randn(2,5)
• Concatenation: join small (compatible)
matrices to make bigger ones:
B = [A A-2; A*2 A/4]

• Deleting rows and columns:
B(:,2) = [ ]
Matrices (con‘t…)

Some useful commands:
zeros(n)             returns a n x n matrix of zeros
zeros(m,n)           returns a m x n matrix of zeros

ones(n)              returns a n x n matrix of ones
ones(m,n)            returns a m x n matrix of ones

size (A)             for a m x n matrix A, returns the row vector [m,n]
containing the number of rows and columns in
matrix.

length(A)            returns the larger of the number of rows or
columns in A.
Matrices (con‘t…)
more commands
Transpose                  B = A‘
Identity Matrix            eye(n)  returns an n x n identity matrix
eye(m,n)  returns an m x n matrix with ones on the main
diagonal and zeros elsewhere.
C=A–B
Scalar Multiplication      B = A, where  is a scalar.
Matrix Multiplication      C = A*B
Matrix Inverse             B = inv(A), A must be a square matrix in this case.
rank (A)  returns the rank of the matrix A.
Matrix Powers              B = A.^2  squares each element in the matrix
C = A * A  computes A*A, and A must be a square matrix.
Determinant               det (A), and A must be a square matrix.
A, B, C are matrices, and m, n,  are scalars.
Solutions to Systems of Linear Equations
•   Example: a system of 3 linear equations with 3 unknowns (x1, x2, x3):
3x1 + 2x2 – x3 = 10
-x1 + 3x2 + 2x3 = 5
x1 – x2 – x3 = -1

Let :
3 2 1               x1         10 
A   1 3 2        x   x2       b5 
                                
 1  1  1
                    x3           1
 
 

Then, the system can be described as:

Ax = b
Solutions to Systems of Linear Equations
(con‘t…)
• Solution by Matrix Inverse:            • Solution by Matrix Division:
Ax = b                                    The solution to the equation
A-1Ax = A-1b                                            Ax = b
x = A-1b                                  can be computed using left division.
• MATLAB:                                  MATLAB:
>> A = [ 3 2 -1; -1 3 2; 1 -1 -1];        >> A = [ 3 2 -1; -1 3 2; 1 -1 -1];
>> b = [ 10; 5; -1];                      >> b = [ 10; 5; -1];
>> x = inv(A)*b                           >> x = A\b
x=                                        x=
-2.0000                                   -2.0000
5.0000                                    5.0000
-6.0000                                   -6.0000
x1 = -2, x2 = 5, x3 = -6                        x1 = -2, x2 = 5, x3 = -6

NOTE:
left division: A\b  b  A           right division: x/y  x  y
Polynomials

• The polynomials are represented by their coefficients in MATLAB.
• Consider the following polynomial:
A(s) = s3 + 3s2 + 3s + 1
• For s is scalar: use scalar operations
– A = s^3 + 3*s^2 + 3*s + 1;
• For s is a vector or a matrix: use array or element by element
operation
– A = s.^3 + 3*s.^2 + 3.*s + 1;
• function polyval(a,s): evaluates a polynomial with coefficients in
vector a for the values in s.
Polynomials (con‘t…)

•   Example:
>> s = linspace (-5, 5, 100);
>> coeff = [ 1 3 3 1];
>> A = polyval (coeff, s);
A(s) = s3 + 3s2 + 3s + 1       >> plot (s, A),
>> xlabel ('s')
>> ylabel ('A(s)')
Polynomials (con‘t…)
Operation         MATLAB command         Description
Addition          c=a+b                  sum of polynomial A and B, the coefficient vectors must be the
same length.
Scalar Multiple   b = 3*a                multiply the polynomial A by 3
Polynomial        c = conv(a,b)          returns the coefficient vector for the polynomial resulting from
Multiplication                           the product of polynomial A and B.
Polynomial        [q,r] = deconv(a,b)    returns the long division of A and B.
Division                                 q is the quotient polynomial coefficient, and r is the remainder
polynomial coefficient.
Derivatives       polyder(a)             returns the coefficients of the derivative of the polynomial A
polyder(a, b)          returns the coefficients of the derivative of the product of
polynomials A and B.
[n,d]=polyder(b,a)     returns the derivative of the polynomial ratio B/A, represented
as N/D
Find Roots        roots(a)               returns the roots of the polynomial A in column vector form.
Find              poly(r)                returns the coefficient vector of the polynomial having roots r.
Polynomials
NOTE: a, b, and c are polynomial coefficient vectors whose corresponding polynomials are A, B, and C
Programming with MATLAB:

• Files that contain code in the MATLAB
language are called M-files. You create M-
files using a text editor, then use them as
you would any other MATLAB functions or
command. There are two types of M-files:
Scripts and Functions.
Scripts
• Scripts: a bunch of code grouped together; doesn‘t accept
argument or return output.
• Example: create a file called magicrank.m that calculates
the rank of magic squares:

r = zeros(1,20);
for n = 1:20
r(n) = rank(magic(n));
end
bar(r)

• Add the file into search path and type the statement:
magicrank
M-Files
So far, we have executed the commands in the command window.
But a more practical way is to create a M-file.

• The M-file is a text file that consists a group of
MATLAB commands.
• MATLAB can open and execute the
commands exactly as if they were entered at
the MATLAB command window.
• To run the M-files, just type the file name in
the command window. (make sure the current
working directory is set correctly)
All MATLAB commands are M-files.
M-files  Functions

• M-files are macros of MATLAB commands that are stored as
ordinary text files with the extension "m", that is filename.m
• example of an M-file that defines a function, create a file in your
working directory named yplusx.m that contains the following
commands:
Write this file:
function z = yplusx(y,x)
z = y + x;
-save the above file(2lines) as yplusx.m
x = 2; y = 3;
z = yplusx(y,x)  5
• Get input by prompting on m-file:
T = input('Input the value of T: ')
Functions:
• Functions are M-files that can accept input
arguments and return output arguments.
The name of the M-file and of the function
should be the same.
• For example, the M-file ―rank.m‖ is
available in the directory:
~toolbox/matlab/matfun, you can see the
file with
type rank
Flow Control:
MATLAB has following flow controls:
•     If statement
•     Switch statement
•     For loops
•     While loops
•     Continue statement
•     Break statement
if … elseif … else … end
If A > B
‗greater‘
elseif A < B
‗less‘
elseif A = = B
‗equal‘
else error(‗Unexpected situation‘)
end
for … end

for i = 1:m
for j = 1:n
H(i,j) = 1/(i+j)
end
end
Suppressing Output:
If you simply type a statement and press Enter,
MATLAB automatically displays the results on
screen. If you end the line with a semicolon ‗;‘,
MATLAB performs the computation but does not
display any result.

• Example: C = randn(5,1) v.s.
C = randn(5,1);
Graphics: plot
x = 0 : .01 : 2*pi;
y = sin(x);
plot(x,y)

y2 = sin(x-.25)
y3 = sin(x-.5)
plot(x,y,x,y2,x,y3)
Plot Continue…
• Adding plots to an existing graph:
hold on

• Multiple plots in one figure:
subplot
Plotting (con‘t…)
•   Example: (polynomial function)
plot the polynomial using linear/linear scale, log/linear scale, linear/log scale, & log/log scale:
y = 2x2 + 7x + 9

% Generate the polynomial:
x = linspace (0, 10, 100);
y = 2*x.^2 + 7*x + 9;

% plotting the polynomial:
figure (1);
subplot (2,2,1), plot (x,y);
title ('Polynomial, linear/linear scale');
ylabel ('y'), grid;
subplot (2,2,2), semilogx (x,y);
title ('Polynomial, log/linear scale');
ylabel ('y'), grid;
subplot (2,2,3), semilogy (x,y);
title ('Polynomial, linear/log scale');
xlabel('x'), ylabel ('y'), grid;
subplot (2,2,4), loglog (x,y);
title ('Polynomial, log/log scale');
xlabel('x'), ylabel ('y'), grid;
Plotting (con‘t…)
Plotting (con‘t…)
•    Adding new curves to the existing graph:
•    Use the hold command to add lines/points to an existing plot.
–      hold on – retain existing axes, add new curves to current axes. Axes are
rescaled when necessary.
–      hold off – release the current figure window for new plots
      Grids and Labels:

Command                        Description
grid on                        Adds dashed grids lines at the tick marks
grid off                       removes grid lines (default)
grid                           toggles grid status (off to on, or on to off)
title (‗text‘)                 labels top of plot with text in quotes
xlabel (‗text‘)                labels horizontal (x) axis with text is quotes
ylabel (‗text‘)                labels vertical (y) axis with text is quotes
text (x,y,‘text‘)              Adds text in quotes to location (x,y) on the current axes, where (x,y) is in
units from the current plot.

color of the point or curve   Marker of the data points   Plot line styles
Symbol    Color              Symbol    Marker             Symbol      Line Style
y       yellow                                             –
.                                     solid line
m      magenta
o                            :        dotted line
c       cyan
x                           –.        dash-dot line
r        red
+         +                  ––        dashed line
g       green
*         
b       blue
w       white                s         □

k       black                d         ◊
v         
^         
h      hexagram
3D - Plotting example
• x=[0:10]; y=[0:10]; z=x‘*y;
• mesh(x,y,z); title(‗3-D Graph‘);
More Plotting
• Old plot got stomped
– To open a new graph, type ‗figure‘
• Multiple data sets:
– Type ‗hold on‘ to add new plot to current
graph
– Type ‗hold off‘ to resume stomping
– title(‗apples over oranges‘)
– xtitle(‗apples‘)
– ytitle(‗oranges‘)
Useful Commands
• Single quote is transpose
• %        same as // comment in C, Java
No /* block comments */ (annoying)
• ;        suppresses printing
• More:
max(x)                   min(x)
mean(x)          median(x)
abs(x)           dot(x,y)
cross(x,y)       flops (flops in this session)
Useful Constants
• Inf     infinity
• NaN            Not and number (div by
zero)
• eps                 machine epsilon
• ans                 most recent
• pi             3.14159….
• i and j        Matlab supports imaginary
numbers!
Final note
• Wrong:
for x = 1:10
for y = 1:10
foo(x,y) = 2 * bar(x,y)
end
end
• Right:
foo = 2 * bar;

• Matlab is optimized for vectorization
Ανθρώπινη όραση, αντίλήψη ιατρικής εικόνας
Ανηίληυη εικόνας: Το ανθρώπινο ζύζηημα

Retina - the ―focal plane array‖ on the
back surface of the eye the detects and
measures light.

Photoreceptors - the nerves in the
retina that detect light.

Fovea - a small region in the retina with
high spatial resolution.

Blind spot - a small region in the retina
where the optic nerve is located that
has no photoreceptors.

Long, medium, and short cones - the three specific types of codes used to
produce color vision. These cones are sensitive to long (L or red) wavelengths,
medium (M or green) wavelengths, and short (S or blue) wavelengths.
Ανηίληυη εικόνας: Βαζικές Καηεσθύνζεις

Goals in the field of perception:

– Understand contrast and how humans detect changes in images.

– Understand photometric properties of the physical world.

– Understand the percept of ―color‖.

– Learn how to use this understanding to design imaging systems.
Ανηίληυη εικόνας: Βαζικά μεγέθη

A Simple Visual Stimulus
A single uniform dot of luminance L in a large uniform background of
luminance LB.
LB

L

How much difference is necessary for a ―standard observer‖ to notice the
difference between L and LB?

– The just noticeable difference (JND) is the difference that allows an observer
to detect the stimulus 50% of the time.
– ΔJND is the difference in L and LB required to achieve a just noticable
difference.
Ανηίληυη εικόνας ζηη Ιαηρική

“Image perception is intimately connected to the psychophysical
properties of the visual system of the observer, i.e. the ability of the
observer to respond to low-contrast and fine-detail stimuli.”

In medical imaging, information is transferred to the observer in two
steps:

1. Data acquisition and image formation in the imaging system
2. Processing and display of image data

Visual performance in medical imaging can be divided into three
categories:

1. Detection of an abnormality (presence)
2. Recognition of an abnormality (shape and size)
3. Identification (relation to likely disease) of the abnormality
Ανηίληυη εικόνας ζηην Ιαηρική
μεηαθέρει με μια μη-
γραμμική διαδικαζία ένα
περιοριζμένο σποζύνολο                          Αληηθείκελο
από ηα αρτικά δεδομένα
(διεθθαρμένο από θόρσβο                                       Ερκελεία/
και «ηετνήμαηα» ), ώζηε να                                    Πιεροθορία
δημιοσργηθεί η ηελική
εικόνα.

Παραηερεηής

Εηθόλα

Σύζηεκα
απεηθόληζες                                             Ο παραηηρηηής
«θιληράρει» ηα δεδομένα
ηης εικόνας (ζηο τώρο
και ζηο τρόνο) και ζηη
ζσνέτεια ερμηνεύει ηο
περιετόμενό ηης με βάζη
ηης οπηικής
πληροθορίας αλλά και
ηης εμπειρίας ηοσ.
Medical Image Analysis
Αληίιευε θαη δηερκελεία ηαηρηθής εηθόλας
Σε ασηή ηελ περίπηφζε ηο «αληηθείκελο» είλαη κηα ηοκή ηοσ
αλζρφπίλοσ ζώκαηος ε οποία έτεη κηα ποιύπιοθε θσζηθο-
τεκηθή δοκή.

Το ζύζηεκα απεηθόληζες κεηαθέρεη κε κηα κε-γρακκηθή
δηαδηθαζία έλα περηορηζκέλο σποζύλοιο από ηα αρτηθά
δεδοκέλα (δηεθζαρκέλο από ζόρσβο θαη «ηετλήκαηα» ), ώζηε
λα δεκηοσργεζεί ε ηειηθή εηθόλα.

Ο παραηερεηής «θηιηράρεη» ηα δεδοκέλα ηες εηθόλας (ζηο τώρο
θαη ζηο τρόλο) θαη ζηε ζσλέτεηα ερκελεύεη ηο περηετόκελό ηες
κε βάζε ηες οπηηθής πιεροθορίας αιιά θαη ηες εκπεηρίας ηοσ.

Ο παραηερεηής αποθαζίδεη θαη παρέτεη πιεροθορίες ζτεηηθά κε
ηο απεηθοληδόκελο αληηθείκελο ...
Ανηίληυη εικόνας ζηην Ιαηρική

•The human brain instinctively searches for geometrical
patterns of objects in an image, such as border lines —
preferring well-known figures such as circles, rectangles and triangles —
and tries to form those objects.

•Τhe viewing distance must be long enough, so that the
pixels cannot be distinguished individually and thus be
seen as squares. Otherwise, the faculty of vision gives priority to
identification of the square-shaped objects rather than to combining them
into a larger object within the image.

•The resolution properties of the eye are dependent on
both the contrast level and the light intensity impinging
on the eye. The resolution drops from ~7line pairs per mm for film-
lightbox combination (x-rays) to less than 3 for a common CRT monitor!!!
Ανηίληυη εικόνας ζηην Ιαηρική : Γνώμες…

"There are times when an experienced
physician experienced physician sees
a visible lesion clearly and times when
he does not. This is the baffling
problem, apparently partly visual and
partly psychologic. They constitute the
still unexplained human equation in
diagnostic procedures.” Henry Garland,
M.D., 1959

“In the last ten years, our basic
knowledge of physics of
increased to such an extent that
it cries out to be linked with
observer performance studies.”
Kurt Rossmann, 1974
Ανηίληυη εικόνας ζηην Ιαηρική : Παραδείγμαηα

Medical Image Perception: Changes in Image display can lead to
errors in interpretation. In this example the “subtle” abnormality
(lung cancer) can be missed in the left image due to poor
display….
Ανηίληυη εικόνας ζηην Ιαηρική : Παραδείγμαηα

Psychophysical and computer automated evaluation and optimization of
image compression

Task: To evaluate the effect of different
image compression algorithms on
visual performance in clinically relevant
tasks of simulated lesions embedded in
real x-ray coronary angiograms.

The evaluation is often done using
human observers (physicians and non-
physicians) reviewing the test-images.
Ανηίληυη εικόνας ζηην Ιαηρική : Παραδείγμαηα
Computer Observer Models
Psychophysical experiments for evaluation of
image acquisition or image processing
technique are time consuming and costly.
For this reason, there has been a growing
computer
interest in developing
observer models that reliably
reproduce human performance
More recent work concentrates on extending
model observers to visual tasks that include
signals that vary in shape and size.
These latter tasks are more representative of
the real clinical scenario where the lesion has
a variety of shapes and sizes
C.K. Abbey and F.O. Bochud, “Modeling visual detection tasks in correlated image noise with linear
model observers,” The Handbook of Medical Imaging: Volume 1, Progress in Medical Physics and
Psychophysics, (Harold Kundel, Ed.), pp. 629-654, 2000.
Ανηίληυη εικόνας ζηην Ιαηρική : Χαρακηηριζμός
ποιόηηηας εικόνας

• The most common method to assess medical image quality is the use
of Receiver Operating Characteristic (ROC) studies.
• In these studies doctors rate their confidence about the presence of a
given abnormaility in a set of patient images.
• The ROC curve plots the hit rate as a function of the false alarm rate.
• From the ROC curve a bias-free measure of performance known as
the area under the curve (AOC) is derived.
• Many ROC studies rely on establishing ―truth‖ (the gold standard)
about lesion absence/presence on the agreement of a panel of experts
(consensus expert committees).

Eckstein, M.P., Wickens, T.D., Aharonov G., Ruan G., Morioka C.A., Whiting, J.S., Quantifying the limitations of the use of
consensus expert committees in ROC studies, Proceedings SPIE Image Perception, 3340, 128-134, (1998)
Medical Image Perception: ROC curves

Sensitivity: Is the probability of a positive test among patients with a
disease. For example if a diagnostic test was positive in 85 out of 100
patients previously diagnosed with breast cancer (in 15 cases the test
was negative), then the sensitivity of the test is 85%.

Specificity: It is the probability of a test being negative among ―healthy‖
patients. Again, the number of negative tests (among patients without the
disease) has to be divided with the total number of tests to calculate the
specificity.
Medical Image Perception: ROC curves

For example, if we want to describe a diagnostic test foe breast cancer
(e.g. a CAD systems that detects malignant masses) and we get the
following results:

Patients with breast cancer   Patient without breast cancer

Positive test        89 (true positives)           40 (false positives)
Negative test        11 (false negatives)          60 (true negatives)

Then we can calculate :

Sensitivity = 89 / 100 = 89%, Specificity = 60 /100= 60%
Medical Image Perception: ROC curves

• If we redefine what is ―positive‖ or ―negative‖ (e.g. change the threshold
at which a detected mass is classified as abnormal), the values for the
specificity and sensitivity will change.
• By doing that for several thresholds we can plot the values of sensitivity
vs. (1-specificity).
• This is called a receiver-operating characteristic curve (ROC).
• Examining different thresholds represent the trade-off between
sensitivity and specificity, or between false positives and negatives.
• The area under the ROC curve is measure of good the test is in
discriminating between ―healthy‖ and non-―healthy‖ patients and can be
used for comparison with other tests.
Medical Image Perception: ROC curves

Patients with breast cancer                Patient without breast cancer
Malignant                   true positives                             false positives
Benign                      false negatives                            true negatives

Βy changing the threshold at which
(e.g. based on texture) a mass was                          100
90
classified as malignant, we got the                          80
following values for the specificity and                     70

sensitivity
sensitivity:                                                 60
50
Threshold that   Sensitivity   Specificity                  40
30
defines
20
malignancy                                                 10                                    ROC

Threshold 1        60%           98%                        0
0   20      40        60    80     100
Threshold 2        80%           94%                                        1-specificity

Threshold 3        90%           85%                      the ROC curve shows that the test is
Threshold 4        95%           75%                      adequate for discriminating between
benign and malignant masses:
Image Analysis and Processing

How to…
Comp. Vision    • Represent
• Process / Prepare
• Handle
• Recognize
• Retrieve
Databases

…images / image objects
Representation

•   Digital Images
•   Color Spaces
•   Gray Images
•   Binary Images
•   Geometrical Properties
Representation

•   Digital Images
•   Color Spaces
•   Gray Images
•   Binary Images
•   Geometrical Properties
Processing

•   Filters
•   Edges
•   Geometric Primitives
• Lines, Circles
Processing

Segmentation
Handling:

•   Image Data Representation
•   Image / Video Formats
• JPEG
• GIF
• MPEG
(this lecture will NOT cover this !)
Low Level Object Handling:

•   Object representation
Object Recognition:
• Color, Texture, Shape
Object Recognition:

• Applications
•   Character recognition
•   Face Recognition
•   Shape Recognition
•   Motion, Movement Detection
•   Behaviour Analysis
•   …
Image Processing Toolbox (2)

• Grayscale Image
(row x col)

• Binary Image
(row x col)
= Grayscale with
2 levels
Image Processing Toolbox (3)

• RGB Image
(row X col X 3)
Image Processing Toolbox (4)
• Indexed Image (2 matrices – colormap
and index)
Image Processing Toolbox (5)
• Reading an Image and storing it in matrix
I:
– [I,map] = imread(‘pout.tif’); For indexed
images
• Deals with many formats
– JPEG, TIFF, GIF, BMP, PNG, PCX, … more can
be added from Mathworks Central File
Exchange
• http://www.mathworks.com/matlabcentral/fileexcha
nge
Image Processing Toolbox (6)
• After an image has been read we can
convert from one type to another:
– ind2gray, gray2ind, rgb2gray, gray2rgb,
rgb2ind, ind2rgb
• Images are read into uint8 data type. To
manipulate the pixel values, they have to
be first converted to double type using the
double(I) function.
Image Processing Toolbox (7)
• Pixel values are accessed as matrix
elements.
– 2D Image with intensity values:            I(row,col)
– 2D RGB images I(row,col,color)
• Color : Red = 1; Green = 2 ; Blue = 3
• Displaying images
– >>figure, imshow(I)
• Displaying pixel position and intensity
information
– pixval on
Image Processing Toolbox (8)
• All arithmetic operations performed on
matrices may be performed on images
• After processing, an image matrix can be
written to an output image file with the
imwrite function
– imwrite(I,map,’filename’,’fmt’)
• Without the map argument, the image data
is supposed to be grayscale or RGB.
• The format ‗fmt‘ needs to support the
particular type of image
Let‘s move on to…

Image Filtering
Image Filtering

• Many basic image processing techniques are based on
convolution.
• In a convolution, a convolution filter is applied to every pixel
to create a filtered image I*(x, y):

      
I * ( x, y )  I ( x , y ) * W ( x , y )      I (u, v)W (u  x, v  y)
u   v  
Image Filtering
• Example: Averaging filter:
y
0     0     0        0      0

0    1/9   1/9       1/9    0

0    1/9   1/9       1/9    0

0    1/9   1/9       1/9    0
x
0     0     0        0      0

      
I * ( x, y )  I ( x , y ) * W ( x , y )       I (u, v)W (u  x, v  y)
u   v  
Image Filtering
• Image:                                       • Filter:

1     6      3     2      9                     1/9       1/9   1/9

2     11     3     10     0                     1/9       1/9   1/9

5     10     6     9      7                     1/9       1/9   1/9

3     1      0     2      8
•New value:
4     4      2     9     10                 •6/9 + 9/9 + 7/9 + 0/9 + 2/9 + 8/9 +
2/9 + 9/9 + 10/9 = 5.889

         
I * ( x, y )  I ( x , y ) * W ( x , y )      I (u, v)W (u  x, v  y)
u   v  
Image Filtering
• Grayscale Image:      • Averaging Filter:

1   6    3   2    9     1/9   1/9   1/9

2   11   3   10   0     1/9   1/9   1/9

5   10   6   9    7     1/9   1/9   1/9

3   1    0   2    8

4   4    2   9    10
Image Filtering
• Original Image:               • Filtered Image:

1     6     3     2   9          0    0    0   0      0
1/9   1/9   1/9

2   11
1/9 1/9
3
1/9
10    0          0    5               0

5   10      6     9   7          0                    0
1/9 1/9     1/9

3     1     0     2   8          0                    0

4     4     2     9   10         0    0    0   0      0

• value = 11/9 + 61/9 + 31/9 + 21/9 + 111/9 + 31/9
+ 51/9 + 101/9 + 61/9 = 47/9 = 5.222
Image Filtering
• Original Image:               • Filtered Image:

1    6      3     2     9        0    0   0    0       0
1/9    1/9   1/9

2    11
1/9
3   10
1/9 1/9
0        0    5   7            0

5    10     6     9     7        0                     0
1/9   1/9   1/9

3    1      0     2     8        0                     0

4    4      2     9     10       0    0   0    0       0

• value = 61/9 + 31/9 + 21/9 + 111/9 + 31/9 + 101/9
+ 101/9 + 61/9 + 91/9 = 60/9 = 6.667
Image Filtering
• Original Image:              • Filtered Image:

1    6   3    2   9            0   0    0   0    0

2   11   3   10   0            0   5    7   5    0

5   10   6    9   7            0   5    6   5    0

3    1   0    2   8            0   4    5   6    0

4    4   2    9   10           0   0    0   0    0

•Now you can see the averaging (smoothing) effect of the 33
filter that we applied.
Image Filtering
• More common: Gaussian Filters

x2  y2
1            
W ( x, y )  G( x, y )              e        2 2
2   2

• implement decreasing influence
by more distant pixels
•1   •4   •7   •4   •1

•4   •16 •26 •16    •4

• Discrete version:                  1/273             •7   •26 •41 •26    •7

•4   •16 •26 •16    •4

•1   •4   •7   •4   •1
Image Filtering
• Effect of Gaussian smoothing:

original         33             99   1515
Different Types of Filters
• Smoothing can reduce noise in the image.
• This can be useful, for example, if you want to find regions
of similar color or texture in an image.

• However, there are different types of noise.
• For so-called ―salt-and-pepper‖ noise, for example, a
median filter can be more effective.
Median Filter
• Use, for example, a 33 filter and move it across the image
like we did before.
• For each position, compute the median of the brightness
values of the nine pixels in question.
– To compute the median, sort the nine values in
ascending order.
– The value in the center of the list (the fifth
value) is the median.
• Use the median as the new value for the center pixel.
Median Filter
• Advantage of the median filter: Capable of eliminating
outliers such as the extreme brightness values in salt-and-
pepper noise.

• Disadvantage: The median filter may change the contours
of objects in the image.
Median Filter

• original image     • 33 median   • 77 median
Image Scaling

This image is too big to
fit on the screen. How
can we reduce it?

How to generate a half-
sized version?
Image sub-sampling

1/8
1/4

Throw away every other row and
column to create a 1/2 size image
Why does this look so crufty?
• Called nearest-neighbor sampling
Even worse for synthetic
images
Sampling and the Nyquist rate

• Aliasing can arise when you sample a
continuous signal or image
– Demo applet
http://www.cs.brown.edu/exploratories/freeSoftware/repository/edu/brown/c
s/exploratories/applets/nyquist/nyquist_limit_java_plugin.html

– occurs when your sampling rate is not high
enough to capture the amount of detail in your
image
– formally, the image contains structure at
different scales
• called ―frequencies‖ in the Fourier domain
Sampling and the Nyquist rate

– the sampling rate must be high enough to
capture the highest frequency in the image
• To avoid aliasing:
– sampling rate > 2 * max frequency in the image
• i.e., need more than two samples per period
– This minimum sampling rate is called the
Nyquist rate
Subsampling with Gaussian
pre-filtering

G 1/8
G 1/4

Gaussian 1/2
• Solution: filter the image, then subsample
– Filter size should double for each ½ size
reduction.
Some times we want many
resolutions

• Known as a Gaussian Pyramid [Burt and Adelson, 1983]
– In computer graphics, a mip map [Williams, 1983]
– A precursor to wavelet transform
• Gaussian Pyramids have all sorts of applications in computer
vision
Gaussian pyramid
construction

• Repeat
– Filter
– Subsample
• Until minimum resolution reached
– can specify desired number of levels (e.g., 3-level
pyramid)
Subsampling with
bilinear pre-filtering

BL 1/8
BL 1/4

Bilinear 1/2
Bilinear interpolation
• A common method for
resampling images
….
Matlab homepage (news & more):
http://www.mathworks.com/

online tutorials:
http://www.engin.umich.edu/group/ctm/
http://www.math.mtu.edu/~msgocken/intro/intro.html

you can find all this at:
http://www.soton.ac.uk/~jowa/teaching.html

DOCUMENT INFO
Shared By:
Categories:
Tags:
Stats:
 views: 169 posted: 11/9/2011 language: English pages: 130
How are you planning on using Docstoc?