Description
AMATH 482-582 Homework 1: An ultrasound problem
You dog fluffy swallowed a marble. The vet suspects that it has now worked its way into the intestines. Using
ultrasound, data is obtained concerning the spatial variations in a small area of the intestines where the marble
is suspected to be. Unfortunately, fluffy keeps moving and the internal fluid movement through the intestines
generates highly noisy data.
Do you want your dog to live? In order to save your dog’s life you must located and compute the trajectory
of the marble. Go to the class webpage and download: Testdata.mat. This containes 20 rows of data for 20
different measurements that were taken in time.
1. Through averaging of the spectrum, determine the frequency signature (center frequency) generated by
the marble.
2. Filter the data around the center frequency determined above in order to denoise the data and determine
the path of the marble. (use plot3 to plot the path once you have it)
3. Where should an intense acoustic wave be focused to breakup the marble at the 20th data measurement.
Good luck, and I hope your dog doesn’t die.
The following code will help you get started in analyzing the data. It also tells you the spatial and spectral
resolution of your ultrasound equipment. (NOTE: the reason for the close all command before isosurface is
that isosurface doesn’t seem to clear the previous imagine before plotting a new one)
clear all; close all; clc;
load Testdata
L=15; % spatial domain
n=64; % Fourier modes
x2=linspace(-L,L,n+1); x=x2(1:n); y=x; z=x;
k=(2*pi/(2*L))*[0:(n/2-1) -n/2:-1]; ks=fftshift(k);
[X,Y,Z]=meshgrid(x,y,z);
[Kx,Ky,Kz]=meshgrid(ks,ks,ks);
for j=1:20
Un(:,:,:)=reshape(Undata(j,:),n,n,n);
close all, isosurface(X,Y,Z,abs(Un),0.4)
axis([-20 20 -20 20 -20 20]), grid on, drawnow
pause(1)
end
AMATH 482-582 Homework 2: G´abor transforms
Part I
In this homework, you will analyze a portion of Handel’s Messiah with time-frequency analysis.
To get started, you can use the following commands (note that Handel is so highly regarded, that
MATLAB has a portion of his music already in MATLAB! I’m sure they’ll add some Beyonce and
T-Swift at some point, but you’ll have to suffer Handel for now):
load handel
v = y’/2;
plot((1:length(v))/Fs,v);
xlabel(’Time [sec]’);
ylabel(’Amplitude’);
title(’Signal of Interest, v(n)’);
This code gives you plots the portion of music you will analyze. To play this back in MATLAB, you
can use the following commands:
p8 = audioplayer(v,Fs);
playblocking(p8);
This homework is rather open ended in the sense that I want you to explore the timefrequency signature of this 9 second piece of classic work.
Things you should think about
doing:
1. Through use of the G´abor filtering we used in class, produce spectrograms of the piece of work.
2. Explore the window width of the G´abor transform and how it effects the spectrogram.
3. Explore the spectrogram and the idea of oversampling (i.e. using very small translations of
the G´abor window) versus potential undersampling (i.e. using very course/large translations
of the G´abor window).
4. Use different G´abor windows. Perhaps you can start with the Gaussian window, and look to
see how the results are effected with the Mexican hat wavelet and a step-function (Shannon)
window.
Don’t be cheap on time here, i.e. don’t be lame. This is an opportunity for you to have a creative
and open ended experience with MATLAB and data analysis. Please do a nice job writing things up
and labeling your resulting plots. I believe this homework can be a really engaging and educational
experience if you devote some time to it.
Part II
From the website, download the two files music1.wav and music2wav. These files play the
song Mary had a little lamb on both the recorder and piano. These are .wav files I generated using
my iPhone. To important and convert them, use the following commands for both pieces (NOTE:
basically both pieces are converted to a vector representing the music, thus you can easily edit the
music by modifying the vector).
Figure 1: Music scale along with frequency of each note in Hz
tr_piano=16; % record time in seconds
y=wavread(’music1’); Fs=length(y)/tr_piano;
plot((1:length(y))/Fs,y);
xlabel(’Time [sec]’); ylabel(’Amplitude’);
title(’Mary had a little lamb (piano)’); drawnow
p8 = audioplayer(y,Fs); playblocking(p8);
figure(2)
tr_rec=14; % record time in seconds
y=wavread(’music2’); Fs=length(y)/tr_rec;
plot((1:length(y))/Fs,y);
xlabel(’Time [sec]’); ylabel(’Amplitude’);
title(’Mary had a little lamb (recorder)’);
p8 = audioplayer(y,Fs); playblocking(p8);
1. Through use of the G´abor filtering we used in class, reproduce the music score for this simple
piece. See Fig. 1 which has the music scale in Hertz. (note: to get a good clean score, you
may want to filter out overtones… see below).
2. What is the difference between a recorder and piano? Can you see the difference in the timefrequency analysis? Note that many people talk about the difference of instruments being
related to the timbre of an instrument. The timbre is related to the overtones generated by
the instrument for a center frequency. Thus if you are playing a note at frequency ω0, an
instrument will generate overtones at 2ω0, 3ω0, · · · and so forth.
AMATH 482-582 Homework 3: PCA
On the webpage are movie files (turned into matlab files) created from three different cameras (videos are from
2011). The experiments are an attempt to illustrate various aspects of the PCA and its practical usefulness
and the effects of noise on the PCA algorithms.
• (test 1) Ideal case: Consider a small displacement of the mass in the z direction and the ensuing
oscillations. In this case, the entire motion is in the z directions with simple harmonic motion being
observed (camN 1.mat where N=1,2,3).
• (test 2) noisy case: Repeat the ideal case experiment, but this time, introduce camera shake into the
video recording. This should make it more difficult to extract the simple harmonic motion. But if the
shake isn’t too bad, the dynamics will still be extracted with the PCA algorithms. (camN 2.mat where
N=1,2,3)
• (test 3) horizontal displacement: In this case, the mass is released off-center so as to produce motion
in the x−y plane as well as the z direction. Thus there is both a pendulum motion and a simple harmonic
oscillations. See what the PCA tells us about the system. (camN 3.mat where N=1,2,3)
• (test 4) horizontal displacement and rotation: In this case, the mass is released off-center and
rotates so as to produce motion in the x−y plane, rotation as well as the z direction. Thus there is both
a pendulum motion and a simple harmonic oscillations. See what the PCA tells us about the system.
(camN 4.mat where N=1,2,3)
Explore the PCA method on this problem and see what you find.
AMATH 482-582 Homework 4: Extended Yale Faces B Database – Eigenfaces & Music Genre Identification
Yale Faces B
Download two data sets (ORIGINAL IMAGE and CROPPED IMAGES)
Your job is to perform an SVD analysis of these data sets. Please start with the cropped images and
perform the following analysis.
1. Do an SVD analysis of the images (where each image is reshaped into a column vector and
each column is a new image).
2. What is the interpretation of the U, Σ and V matrices?
3. What does the singular value spectrum look like and how many modes are necessary for good
image reconstructions? (i.e. what is the rank r of the face space?)
4. compare the difference between the cropped (and aligned) versus uncropped images.
This is an exploratory homework. So play around with the data and make sure to plot the different
things like the modes and singular value spectrum. Good luck, and have fun.
Music Classification
Music genres are instantly recognizable to us, wether it be jazz, classical, blues, rap, rock, etc. One
can always ask how the brain classifies such information and how it makes a decision based upon
hearing a new piece of music. The objective of this homework is to attempt to write a code that can
classify a given piece of music by sampling a 5 second clip.
As an example, consider Fig. 1. Four classic pieces of music are demonstrated spanning genres
of rap, jazz, classic rock and classical. Specifically, a 3-second sample is given of Dr. Dre’s Nuthin’
but a ’G’ thang (The Chronic), John Coltrane’s A Love Supreme (A Love Supreme), Led Zeppelin’s
Over The Hills and Far Away (Houses of the Holy), and Mozart’s Kyrie (Requiem). Each has a
different signature, thus begging the question wether a computer could distinguish between genres
based upon such a characterization of the music.
• (test 1) Band Classification: Consider three different bands of your choosing and of different
genres. For instance, one could pick Michael Jackson, Soundgarden, and Beethoven. By taking
5-second clips from a variety of each of their music, i.e. building training sets, see if you can
build a statistical testing algorithm capable of accurately identifying ”new” 5-second clips of
music from the three chosen bands.
• (test 2) The Case for Seattle: Repeat the above experiment, but with three bands from
within the same genre. This makes the testing and separation much more challenging. For
instance, one could focus on the late 90s Seattle grunge bands: Soundgarden, Alice in Chains,
and Pearl Jam. What is your accuracy in correctly classifying a 5-second sound clip? Compare
this with the first experiment with bands of different genres.
• (test 3) Genre Classification: One could also use the above algorithms to simplify broadly
classify songs as jazz, rock, classical etc. In this case, the training sets should be various bands
2 2.5 3 3.5 4 4.5 5
−1
0
1
2 2.5 3 3.5 4 4.5 5
−0.5
0
0.5
2 2.5 3 3.5 4 4.5 5
−0.2
0
0.2
2 2.5 3 3.5 4 4.5 5
−0.2
0
0.2
time (seconds)
Figure 1: Instantly recognizable, these four pieces of music are (in order of top to bottom): Dr.
Dre’s Nuthin’ but a ’G’ thang (The Chronic), John Coltrane’s A Love Supreme (A Love Supreme),
Led Zeppelin’s Over The Hills and Far Away (Houses of the Holy), and Mozart’s Kyrie (Requiem).
Illustrated is a 3-second clip from time 2 seconds to 5 seconds of each of these songs.
within each genre. For instance, classic rock bands could be classified using sounds clips from
Zep, AC/DC, Floyd, etc. while classical could be classified using Mozart, Beethoven, Bach,
etc. Perhaps you can limit your results to three genres, for instance, rock, jazz, classical.
WARNING and NOTES: You will probably want to SVD the spectrogram of songs versus the songs
themselves. Interestingly, this will give you the dominant spectrogram modes associated with a given
band. Moreover, you may want to re-sample your data (i.e. take every other point) in order to keep
the data sizes more manageable. Regardless, you will need lots of processing time.
AMATH 482-582 Homework 5: Background Subtraction in Video Streams
Use the Dynamic Mode Decomposition method to take a video clip containing a foreground and
background object and separate the video stream to both the foreground video and a background.
Make several test videos to try the algorithm on. You can use your smart phone to generate appropriate videos.
The DMD spectrum of frequencies can be used to subtract background modes. Specifically,
assume that ωp, where p ∈ {1, 2, . . . , `}, satisfies ∥ωp∥ ≈ 0, and that ∥ωj∥ ∀ j ≠ p is bounded away
from zero.
Thus,
XDMD = bpϕp
e
ωpt
´ ¹¹¹¹¹¹¹¹¸ ¹¹¹¹¹¹¹¹¶
Background Video
+ ∑
j≠p
bjϕj
e
ωj t
´ ¹¹¹¹¹¹¹¹¹¹¹¹¹¸ ¹¹¹¹¹¹¹¹¹¹¹¹¹¶
Foreground Video
(1)
Assuming that X ∈ R
n×m, then a proper DMD reconstruction should also produce XDMD ∈ R
n×m.
However, each term of the DMD reconstruction is complex: bjϕj
exp (ωjt) ∈ C
n×m ∀j, though they
sum to a real-valued matrix.
This poses a problem when separating the DMD terms into approximate
low-rank and sparse reconstructions because real-valued outputs are desired and knowing how to
handle the complex elements can make a significant difference in the accuracy of the results. Consider
calculating the DMD’s approximate low-rank reconstruction according to
XLow-Rank
DMD = bpϕp
e
ωpt
.
Since it should be true that
X = XLow-Rank
DMD + X
Sparse
DMD ,
then the DMD’s approximate sparse reconstruction,
X
Sparse
DMD = ∑
j≠p
bjϕj
e
ωj t
,
can be calculated with real-valued elements only as follows. . .
X
Sparse
DMD = X − ∣XLow-Rank
DMD ∣,
where ∣ ⋅ ∣ yields the modulus of each element within the matrix. However, this may result in X
Sparse
DMD
having negative values in some of its elements, which would not make sense in terms of having
negative pixel intensities.
These residual negative values can be put into a n × m matrix R and then
be added back into XLow-Rank
DMD as follows:
XLow-Rank
DMD ← R + ∣XLow-Rank
DMD ∣
X
Sparse
DMD ← X
Sparse
DMD − R
This way the magnitudes of the complex values from the DMD reconstruction are accounted for,
while maintaining the important constraints that
X = XLow-Rank
DMD + X
Sparse
DMD ,
so that none of the pixel intensities are below zero, and ensuring that the approximate low-rank and
sparse DMD reconstructions are real-valued. This method seems to work well empirically.
NOTE: it is pretty easy to produce a video on your smart phone.