Introduction To The Frequency Domain Alongside R
This article is somewhat experimental on my part. If you lot are unfamiliar amongst frequency domain concepts, my feeling is the best thought would live on to experiment amongst the R commands that I used to generate these screenshots. Understanding why this matters may need to hold off until I write follow upward articles explaining how frequency domain concepts influence fourth dimension serial analysis.
What Is The Frequency Domain?
The basic concept of a frequency domain representation of a signal at nowadays appears to live on good known every bit a number of advances inwards scientific discipline together with electronics. The difficulty in all probability lies inwards the mathematics some the concept. But it is possible to prepare an intuitive agreement of the frequency domain, without recourse to the underlying mathematical formulae.
 |
- A Fourier Transform converts a fourth dimension domain serial into a frequency domain representation.
- An Inverse Fourier Transform converts dorsum a frequency domain representation into a fourth dimension series.
Please Federal Reserve notation that at that topographic point are a number of Fourier transforms, for both continuous fourth dimension together with discrete time, every bit good every bit related operators (Laplace together with z-transform). Within this article, I volition work the Fast Fourier Transform, the "FFT", which is heavily used inwards practice.
The frequency domain representation of a fourth dimension serial is generated yesteryear breaking the serial upward into a laid of underlying sinusoidal serial of a number of dissimilar frequencies. (A sinusoid is a fourth dimension serial that is generated yesteryear a business office that is a combination of sine and/or cosine functions.)
As long every bit nosotros do non enquire awkward questions every bit to how the FFT works, nosotros tin care for it every bit a dark box together with only examine the output empirically. Since I believe that frequency domain concepts are mainly useful for agreement the telephone commutation limitations of what analysis tin live on done on economical fourth dimension series, together with less useful for computational purposes, I do non consider the need to delve likewise deep into the mathematical formulae.
Using The FFT In R
The start thing to do is to generate a sine moving ridge inwards the fourth dimension domain. In the R code above, I create a fourth dimension vector that runs from t=0 to t=63, which is 64 points. The preferred indexing to piece of work amongst is to receive got the start index every bit existence 0, together with non 1. But the R linguistic communication defines the start signal inwards a vector every bit having an index of 1, which makes things messier.
The line:
s <- sin(2*pi*t/16) # Frequency = 1/16generates the entire fourth dimension serial (denoted s) inwards ane occupation of code; this is because the R linguistic communication is vectorised. In some languages, it would live on necessary to run a loop to generate all the values of the fourth dimension series. Also Federal Reserve notation that R uses the "s <-" to assign a value to the variable s. As the plot shows, the fourth dimension serial s repeats every sixteen fourth dimension steps. (Note that I receive got truncated the x-axis yesteryear the work of a xlim option.)
Please Federal Reserve notation that the frequency of this fourth dimension serial is 1/16 (0.0625), which is expressed relative to the fourth dimension steps inwards discrete time. If the fourth dimension steps were ane minute apart, this would check to a frequency of 0.0625 Hertz (1 Hertz = 1/second). However, since nosotros do non know how much fourth dimension each fourth dimension mensuration represents, the frequencies I quote hither receive got no dimension. To set this into context for economics, an annual frequency when nosotros sample information quarterly corresponds to a (dimensionless) frequency of 1/4 (= 0.25), but if nosotros sample information monthly, an annual frequency is 1/12 (=.0833). In other words, the same existent the world frequency maps to dissimilar dimensionless frequencies depending upon the charge per unit of measurement at which you lot sample data.
This fourth dimension serial s volition live on used repeatedly throughout this article.
The screenshot higher upward shows the calculation of the frequency domain representation of s, which is assigned to f_s. When I impress out the start half dozen entries of f_s, we teach a pocket-sized surprise: the values of the vector are at nowadays complex numbers. That is, they receive got an imaginary component, denoted amongst i, where i is the foursquare root of -1. (In electrical applied scientific discipline texts, j is used to stand upward for this number, every bit i was already reserved to refer to electrical flow inwards a circuit.)
Since it is hard to visualise imaginary numbers, I plotted the absolute value of each of the vector components. (The absolute value is to a greater extent than formally the modulus of the complex number.) We tin consider the magnitude of each entry, but nosotros lose the information almost the breakdown of the existent together with imaginary components.
As the nautical chart shows, the magnitude of most entries is nearly goose egg (with non-significant calculation errors), except for the fifth together with 61st entry of the vector, which are -32i and +32i respectively.
Why those entries? If nosotros locomote to a zero-based indexing, those are entries iv together with lx (out of 64). The entry iv is equal to the size of the sample (64) divided yesteryear the frequency (16). The other entry is a complement of that value. Since I did non desire to delve into the mathematics, all I tin state is that this complement has to live on there. I volition demo below what happens if that higher frequency ingredient is deleted.
Also Federal Reserve notation that at that topographic point is an asymmetry. The start entry (index 0) corresponds to a constant signal (or straight off electrical flow (DC) inwards electrical engineering). If nosotros accept the FFT of fourth dimension serial that is 1 everywhere, the FFT is goose egg everywhere except the start entry.
> x = rep(1,64)
> x
[1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
[37] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
> fft(x)
[1] 64+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i
[13] 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i
[25] 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i
[37] 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i
[49] 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i 0+0i
[61] 0+0i 0+0i 0+0i 0+0i
In summary, the frequency domain representation volition receive got a start entry that corresponds to a constant offset, together with and then the remaining entries check to time-varying sinusoids.
Phase Difference
In the code above, I create some other fourth dimension series, called c, which is generated yesteryear the cosine (cos) function. As tin live on seen inwards the plot, the cosine is only similar the sine wave, except that it is leading yesteryear iv periods. This corresponds to a number from trigonometry which most people receive got in all probability forgotten (including the author, who had to validate the result):
cos(x) = sin(x + π/2).
When nosotros compare the frequency domain representations, they receive got the same magnitude at all points. (If nosotros subtract the absolute values of the frequency domain vectors, the deviation is zero.)
The deviation betwixt the 2 frequency domain representations is the breakdown of the real/imaginary parts. For the cosine, the fifth entry is 32, compared to the sine which is -32i.
In general:
- the magnitude (modulus) of a frequency domain entry gives the aAmplitude of the sine wave;
- the relative sizes of the real/imaginary parts of the frequency domain entry determines the phase shift (lead/lag) of the sine wave.
Filtering Influenza A virus subtype H5N1 Noisy Signal
We tin at nowadays facial expression what happens when nosotros brand the fourth dimension serial to a greater extent than complicated. This was done yesteryear adding "noise" which is only some other sinusoid of higher frequency.
The code higher upward creates a racket signal n which at nowadays repeats every 8 periods (frequency = 1/8, which is larger than 1/16). The serial "signal" is the amount of the racket together with the master sine moving ridge (signal = s + n; note that nosotros tin work "=" for assignment every bit good inwards R). The resulting fourth dimension domain signal is shown above.
As ane mightiness expect, the FFT of "signal" at nowadays has 2 peaks, at (zero-based) indices iv together with 8; which check to 64/16, together with 64/8 respectively.
We tin together with then apply an ideal low-pass filter , which only eliminates high frequency components. I only goose egg out the middle of the frequency domain representation of the signal.
We tin together with then convert this filtered signal dorsum to the fourth dimension domain, using the fft(f_filtered,inverse=TRUE) operation. Given a quirk inwards the Definition of the FFT, inwards companionship to select dorsum the truthful fourth dimension domain representation, nosotros receive got to split upward yesteryear the length of the sample menses (in this illustration 64). This agency that the "inverse FFT" is non the "true" Inverse Fourier Transform. The argue for this odd scaling is that the touchstone Definition of the FFT is computationally efficient. (It is called the Fast Fourier Transform for a reason; it was a nifty bound forrard for the efficiency of digital signal processing.)
Therefore, this ideal low-pass filter does what you lot intuitively expect: it wipes out the high frequency racket inwards the signal.
What If We Break Complementarity?
The screenshot higher upward shows what happens if nosotros suspension the complementarity of the frequency domain representation. In it, I take away all of the vector beyond the sixth index. This eliminates the noise, every bit good every bit the top halt of the master lower frequency (frequency = 1/16) signal.
The fourth dimension domain serial subsequently the inverse FFT at nowadays has imaginary parts. If nosotros plot the existent part, nosotros consider that it is a sine wave, but amongst one-half of the amplitude. The modulus is is constant, together with is 0.5, which is one-half of the master amplitude.
What this agency is that the master signal is the amount of 2 of complex sine waves, amongst the imaginary parts cancelling out.
More technically, a FFT on a signal amongst north fourth dimension points is converted into a frequency domain representation sampled at north point. The sample (dimensionless) frequencies are given by: 0, 1/N, 2/N,... (N-1)/N. (This converges towards a frequency of 1, which is a signal that repeats every fourth dimension sample.) For a discrete fourth dimension signal, what happens is that a (dimensionless) frequency that is 1 or higher becomes indistinguishable from a frequency that is inwards the interval [0,1) (numbers that are greater than or equal to 0, but strictly less than 1). For example, if a sinusoid has frequency of 1, it volition repeat the same value each period, together with thus it is indistinguishable from a constant (frequency = 0). (Please Federal Reserve notation that this type of aliasing does non give for continuous fourth dimension signals, which is why this concept is non intuitive, since it is non commonly encountered.) Please Federal Reserve notation that the commons convention is to announce this interval every bit existence [0,2Ï€), but I receive got embedded the 2Ï€ within my Definition of a frequency.
It takes a lilliputian flake of trigonometry, but nosotros volition observe that only the frequencies inwards the hit [0,1/2) tin live on represented inwards the fourth dimension domain; the frequency domain components inwards the hit (1/2,1) are only the complements of the lower frequency signal.
In summary: when nosotros facial expression at the FFT of a fourth dimension domain signal, only the start one-half of the vector genuinely contains information.The minute one-half of the vector only ensures that the inverse functioning returns a real-valued fourth dimension signal.
Appendix: R Code
If you lot re-create this code into a file, the "readline()" command volition halt execution until you lot striking "return". The R programming linguistic communication is available hither (free).#######################
# First plot
# Create 64 menses fourth dimension vector, sine wave
t <- 0:63
# Operation is vectorised, tin do inwards ane line
s <- sin(2*pi*t/16) # Frequency = 1/16
print(s[1:17])
plot(t,s,xlim=c(0,36))
title(main="Sine Wave, frequency=1/16")
abline(h=0,v=c(0,16,32)) # Add mark lines
readline()
#######################
# Second plot
idx <- 0:63
f_s <- fft(s)
# Display elements 1-6: imaginary numbers!
print(f_s[1:6],digits=3)
# Since imaginary, plot the absolute value
plot(idx,abs(f_s))
title("FFT of sine wave")
# Second peak at chemical element #61 (=60 zero-based)
print(f_s[60:62],digits=3)
readline()
#######################
# Third plot
# Create c(t) = moving ridge created yesteryear cosine
c <- cos(2*pi*t/16) #cosine,Frequency = 1/16
# Plot the two
matplot(t,cbind(s,c),type="b",pch=1:2,col=1:2)
legend(x='topright',pch=1:2,legend=c("sine","cosine"),col=1:2)
f_c <- fft(c)
# Display chemical element v of both
mat <- cbind(f_s,f_c)
print(mat[5,],digits=3)
# Note that f_s has fifth chemical element = -32i,
# f_c has fifth chemical element = 32
# However, magnitudes are (almost) identical
print(max(abs(f_c)-abs(f_s)))
readline()
#####################################
# Plot 4: signal & noise
# add together a signal amongst frequncy of 1/8 to the sine
n = sin(2*pi*t/8) # freq=1/8
signal = s + n
plot(t,signal,type="b")
title("Time serial of noisy signal")
readline()
######################################
# Plot 5: FFT of signal & noise:
f_signal = fft(signal)
plot(idx,abs(f_signal))
title("FFT Of Signal & Noise")
readline()
#####################################
# Plot 6:
# Ideal depression overstep filter take away everything beyond
# the sixth element. Note at that topographic point is an index asymmetry.
f_filtered <- f_signal
f_filtered[1+6:(64-6)] = 0
# Need to normalise yesteryear 1/length
filtered <- fft(f_filtered,inverse=TRUE)/64
# Noted: at that topographic point are tiny (1e-16) imaginary parts
# Use Re() to take away them
plot(t,Re(filtered),type="b")
title("Time serial subsequently low-pass filter")
readline()
#####################################
# Plot 7:
# Broken ideal "low-pass" filter
f_filtered <- f_signal
f_filtered[(1+6):64] = 0 # Blow out all of high end
# Need to normalise yesteryear 1/length
filtered <- fft(f_filtered,inverse=TRUE)/64
# Imaginary real-time data!
print(filtered[1:5])
matplot(t,cbind(Re(filtered),abs(filtered)),type="b",pch=1:2,col=1:2)
legend(x='right',pch=1:2,legend=c("Real Part","Modulus"),col=1:2)
title("Real Part of Time serial After Broken Filtering")
(c) Brian Romanchuk 2015
No comments