Convolution matrix reversed impulse response impulse response - - PDF document

Mathematical Tools for Neural and Cognitive Science

Section 3: Linear Shift-invariant Systems

Fall semester, 2018

1

Linear shift-invariant (LSI) systems

• Linearity (previously discussed):

“linear combination in, linear combination out”

• Shift-invariance (new property):

“shifted vector in, shifted vector out”

• Note: These two properties are independent

(think of some examples that have both, one, or neither)

2

LSI system

v

Input

v1 x v4 x v3 x v2 x

L

Output

v1 x v4 x v3 x v2 x + + + + + +

As before, express input as a sum of “impulses”, weighted by elements of x

3

LSI system

v

Input

v1 x v4 x v3 x v2 x

L

Output

v1 x v4 x v3 x v2 x + + + + + +

• Shift-invariance => responses to

impulses are shifted copies of each other

• Linearity => response to x is sum of

responses to impulses, weighted by elements of x

4

LSI system

v

Input

v1 x v4 x v3 x v2 x

L

Output

v1 x v4 x v3 x v2 x + + + + + +

LSI systems are characterized by their “impulse response”

5

Convolution

• Sliding dot product
• Structured matrix
• Boundaries? zero-padding, reflection, circular
• Examples: impulse, delay, average, difference

In Out

+

r

1

r2 r3

+

r0 r1 r2

y(n) = X

k

r(n − k)x(k) = X

k

r(k)x(n − k)

6

Convolution matrix

impulse response reversed impulse response boundaries?

7

In Out

+

r

1

r2 r3

Feedback LSI system

y(n) = X

k

f(n − k)x(k) + X

k

g(n − k)y(k) (For this class, we’ll stick to feedforward (FIR) systems)

• Response depends on input, and

previous outputs

• Infinite impulse response (IIR)
• Recursive => possibly unstable

8

2D convolution

[figure c/o Castleman]

“sliding window”

9

Discrete Sinusoids

More generally: “amplitude” “phase” (radians) “frequency” (radians/sample) “frequency” (cycles/vectorLength) , ω = 2πk/N cos(ωn)

10 20 30 −1 1 10 20 30 −2 −1 1 2

example : A = 1.5, φ = 8π/32

10

Shifting Sinusoids

... via a well-known trigonometric identity: cos(a − b) = cos(a) cos(b) + sin(a) sin(b) We’ll also need conversions between polar and rectangular coordinates: A = p x2 + y2, φ = tan−1(y/x) x = A cos(φ), y = A sin(φ)

x y A φ

A cos(ωn − φ) = A cos(φ) cos(ωn) + A sin(φ) sin(ωn)

11

Shifting Sinusoids

A cos(ωn − φ) = A cos(φ) cos(ωn) + A sin(φ) sin(ωn) Any scaled and shifted sinusoidal vector can be written as a weighted sum of two fixed {sin, cos} vectors!

A sin φ

φ

A cos φ

A

scale factors: fixed cos/sin vectors:

10 20 30 −1 1 10 20 30 −1 1

12

10 20 30 −1 1

A = 1.6, φ = 2π0/12

A cos(ωn − φ) = A cos(φ) cos(ωn) + A sin(φ) sin(ωn)

fixed cos/sin vectors:

10 20 30 −1 1

A = 1.6, φ = 2π1/12

Shifting Sinusoids

10 20 30 −1 1 10 20 30 −1 1 10 20 30 −1

Any scaled and shifted sinusoidal vector can be written as a weighted sum of two fixed {sin, cos} vectors!

13

10 20 30 −1 1

A = 1.6, φ = 2π2/12 A = 1.6, φ = 2π3/12

10 20 30 −1 1

A = 1.6, φ = 2π4/12

10 20 30 −1 1

A = 1.6, φ = 2π5/12

10 20 30 −1 1

A = 1.6, φ = 2π6/12

10 20 30 −1 1

A cos(ωn − φ) = A cos(φ) cos(ωn) + A sin(φ) sin(ωn)

fixed cos/sin vectors:

Shifting Sinusoids

10 20 30 −1 1 10 20 30 −1 1 10 20 30 −1

Any scaled and shifted sinusoidal vector can be written as a weighted sum of two fixed {sin, cos} vectors!

14

x(n) = cos(ωn)

LSI response to sinusoids

(input)

15

x(n) = cos(ωn)

(convolution formula)

LSI response to sinusoids L

16

x(n) = cos(ωn)

inner product of impulse response with cos/sin, respectively (trig identity)

LSI response to sinusoids L

17

x(n) = cos(ωn)

L LSI response to sinusoids

18

x(n) = cos(ωn)

A sin φ

φ

A cos φ

A Rc(ω) Rs(ω)

cr(ω) sr(ω)

(rectangular -> polar coordinates)

LSI response to sinusoids

19

x(n) = cos(ωn)

L

“Sinusoid in, sinusoid out” (with modified amplitude & phase)

LSI response to sinusoids

(trig identity, in the opposite direction)

20

phases add amplitudes multiply

L

“Sinusoid in, sinusoid out” (with modified amplitude & phase) More generally, if input has amplitude and phase , Ax φx then linearity and shift-invariance tell us that

LSI response to sinusoids

21

• Frequency multiples of radians/sample,


(specifically, )

• Construct an orthogonal matrix of sin/cos pairs,

covering different numbers of cycles

The Discrete Fourier transform (DFT)

[details on board...]

• When we apply this matrix to an input vector, think
• f output as paired coordinates
• Common to plot these pairs as amplitude/phase
• For , only need the cosine part

(thus, cosines, and sines)

N/2 − 1 N/2 + 1

22

k=0 k=1 k=2 k=3

F =

k=N/2

cos ( 2πk N n) sin ( 2πk N n)

(plotted sinusoids are continuous, N=32)

Fourier Transform matrix

23

The Fourier family

(we are here)

signal domain frequency domain

The “fast Fourier transform” (FFT) is a computationally efficient implementation of the DFT, requiring Nlog(N) operations, compared to the N2 operations that would be needed for matrix multiplication.

24

x(n) = cos(ωn) NOTE: These dot products are the Fourier transform

• f the impulse response, r(m)!

LSI response to sinusoids

25

⇥ x

L

Fourier & LSI

26

⇥ x

L

Fourier & LSI

note: only 3 (of many) frequency components shown

27

⇥ x

L

Fourier & LSI

note: only 3 (of many) frequency components shown

28

v

Input

v1 x v4 x v3 x v2 x

L

Output

v1 x v4 x v3 x v2 x + + + + + +

LSI systems are characterized by their frequency response, specified by the Fourier Transform of their impulse response ⇥ x

L

Fourier & LSI

29

eiθ = cos(θ) + i sin(θ) eiωn = cos(ωn) + i sin(ωn)

n n

real part: imaginary part:

[on board: reminders of addition/multiplication of complex numbers]

Complex exponentials: “bundling” sine and cosine

30

eiωn

L Complex exponentials: “bundling” sine and cosine

F.T. of impulse response!

31

eiωn

L

F.T. of impulse response!

L

Note: the complex exponentials are eigenvectors!

Complex exponentials: “bundling” sine and cosine

32

convolve with

The “convolution theorem”

33

convolve with Fourier Transform inverse Fourier Transform pointwise multiply by

The “convolution theorem”

34

The “convolution theorem”

convolve with Fourier Transform inverse Fourier Transform pointwise multiply by L FT F

(diagonal matrix)

⇒ F T~ y = ˜ RF T~ x 35

Recap…

• Linear system
• defined by superposition
• characterized by a matrix
• Linear Shift-Invariant (LSI) system
• defined by superposition and shift-invariance
• characterized by a vector (the impulse response)
• OR, characterized by frequency response.

Specifically, the Fourier Transform of the impulse response specifies an amplitude multiplier and a phase shift for each frequency.

36

Discrete Fourier transform (with complex numbers)

where ωk = 2πk N (inverse)

k

rn = 1 N

N−1

X

k=0

˜ rk eiωkn

[on board: why minus sign? why 1/N?]

Redraw with spiral included!

37

Visualizing the (Discrete) Fourier Transform

• Two conventional choices for frequency axis:
• Plot frequencies from k=0 to k=N/2

(in matlab: 1 to N/2-1)

• Plot frequencies from k=-N/2 to N/2-1

(in matlab: use fftshift)

• Typically, plot amplitude (and possibly phase,
• n a separate graph), instead of the real/

imaginary (cosine/sine) components

38

More examples

• constant
• sinusoid (see next slide)
• impulse
• Gaussian - “lowpass”
• DoG (difference of 2 Gaussians) - “bandpass”
• Gabor (Gaussian windowed sinusoid) - “bandpass”

[on board]

39

Example for k=2, N=32 (note indexing and amplitudes):

10 20 30 −1 1 10 20 30 −10 10 10 20 30 −1 1 −10 10 −10 10 k −10 10 −10 10 k 10 20 30 −10 10 (real part) (imag part)

=>

40

What do we do with Fourier Transforms?

• Represent/analyze periodic signals
• Analyze/design LSI systems. In particular,

how do you identify the nullspace?

41

Retinal ganglion cells (1D)

Enroth-Cugell and Robson (1984)

42

Sampling causes “aliasing”

“Aliasing” - one frequency masquerades as another [on board] Sampling process is linear, but many-to-one (non-invertible) Given the samples, it is common/natural to assume, or enforce, that they arose from the lowest compatible frequency...

43

Effect of sampling on the Fourier Transform: Sum of shifted copies |X(ω)| |Xs(ω)|

44

Real-world aliasing

downsample by 2

“Moiré pattern”

45

Pre-filtering to avoid spectral overlap (“aliasing”)

L(ω) L(ω)

lowpass filter, cutoff at π/∆

|X(ω)| |Xs(ω)|

46

Real-world aliasing

, with pre-filtering downsample by 2

47