Python SciPy IIR Filter: Digital Signal Processing

by

Recently, I was working on a project where I needed to process audio signals and remove unwanted noise. As I explored different solutions, I found that SciPy’s IIR filters were incredibly powerful for this task.

In this article, I’ll share my hands-on experience with SciPy’s IIR filters, showing you how to implement them effectively in your Python projects. Whether you’re working with audio processing, sensor data, or any signal processing application, these filters can be game-changers.

Let’s get in!

What is an IIR Filter?

IIR stands for Infinite Impulse Response, a type of digital filter that uses feedback loops. Unlike FIR (Finite Impulse Response) filters, IIR filters can produce an output that theoretically continues indefinitely.

The key characteristics of IIR filters include:

  • Higher efficiency with fewer coefficients
  • Feedback from previous outputs
  • Potential for instability if not designed properly
  • Nonlinear phase response

Read Python Scipy Convolve 2d

Set Up Your Environment

Before we begin implementing IIR filters, make sure you have the necessary libraries installed:

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

If you don’t have SciPy installed yet, you can install it using pip:

pip install scipy numpy matplotlib

Method 1 – Use signal.iirfilter()

The simplest way to create an IIR filter in SciPy is by using the iirfilter() function. This function is versatile and allows you to design various IIR filters based on your requirements.

# Create a sample signal with noise
sample_rate = 1000  # Hz
t = np.linspace(0, 1, sample_rate, endpoint=False)

# Create a 50 Hz signal with high-frequency noise
clean_signal = np.sin(2 * np.pi * 50 * t)
noise = 0.5 * np.sin(2 * np.pi * 250 * t)
signal_with_noise = clean_signal + noise

# Design a 5th order Butterworth lowpass filter with cutoff at 100 Hz
b, a = signal.iirfilter(5, 100/(sample_rate/2), btype='lowpass', 
                        ftype='butter', output='ba')

# Apply the filter
filtered_signal = signal.filtfilt(b, a, signal_with_noise)

# Plot the results
plt.figure(figsize=(10, 6))
plt.plot(t, signal_with_noise, 'b-', alpha=0.5, label='Noisy Signal')
plt.plot(t, filtered_signal, 'r-', label='Filtered Signal')
plt.plot(t, clean_signal, 'g--', alpha=0.7, label='Original Signal')
plt.legend()
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')
plt.title('IIR Filter Example with iirfilter()')
plt.grid(True)
plt.tight_layout()
plt.show()

You can refer to the screenshot below to see the output:

scipy.signal.cheby2 python digital filter

In this example, I’ve designed a 5th-order Butterworth lowpass filter with a cutoff frequency of 100 Hz. The filtfilt() function applies the filter forward and backward, giving zero phase distortion.

Method 2 – Use Specific Filter Design Functions

SciPy provides dedicated functions for common filter types. These are often more intuitive if you know exactly what filter you need.

Butterworth Filter

Butterworth filters have a maximally flat frequency response in the passband.

# Create a sample signal with noise
fs = 1000  # Sampling frequency (Hz)
t = np.linspace(0, 1, fs, endpoint=False)

# Original signal: 40 Hz sine wave + 300 Hz noise
signal_main = np.sin(2 * np.pi * 40 * t)
noise = 0.6 * np.sin(2 * np.pi * 300 * t)
noisy_signal = signal_main + noise

# Design a Butterworth lowpass filter using butter()
cutoff = 100  # Cutoff frequency (Hz)
order = 6
sos = signal.butter(order, cutoff, btype='low', fs=fs, output='sos')

# Apply the filter using sosfiltfilt (zero-phase filtering)
filtered = signal.sosfiltfilt(sos, noisy_signal)

# Plot the original and filtered signals
plt.figure(figsize=(10, 5))
plt.plot(t, noisy_signal, label='Noisy Signal', alpha=0.5)
plt.plot(t, filtered, label='Filtered Signal (Butterworth)', linewidth=2)
plt.plot(t, signal_main, '--', label='Original Signal (40 Hz)', alpha=0.7)
plt.legend()
plt.title('Butterworth Lowpass Filter using signal.butter()')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')
plt.grid(True)
plt.tight_layout()
plt.show()

You can refer to the screenshot below to see the output:

scipy iir filter

I prefer using the SOS (second-order sections) output format as it provides better numerical stability, especially for higher-order filters.

Chebyshev Type I Filter

Chebyshev Type I filters allow ripple in the passband but have a steeper rolloff.

# Design a Chebyshev Type I lowpass filter
b, a = signal.cheby1(4, 1, 100, btype='lowpass', fs=sample_rate)
filtered_signal = signal.lfilter(b, a, signal_with_noise)

Elliptic Filter

Elliptic filters offer the steepest transition between passband and stopband, with ripples in both.

# Design an elliptic bandstop filter to remove 250 Hz noise
sos = signal.ellip(6, 1, 40, [240, 260], btype='bandstop', 
                   fs=sample_rate, output='sos')
filtered_signal = signal.sosfilt(sos, signal_with_noise)

Check out Python Scipy Leastsq

Method 3 – Real-Time Processing with lfilter()

For real-time applications like audio streaming, lfilter() is more appropriate than filtfilt():

# Create a sample signal with high-frequency noise
sample_rate = 1000  # Hz
t = np.linspace(0, 1, sample_rate, endpoint=False)

# Clean 50 Hz signal with 250 Hz noise
clean_signal = np.sin(2 * np.pi * 50 * t)
noise = 0.5 * np.sin(2 * np.pi * 250 * t)
signal_with_noise = clean_signal + noise

# Simulate real-time processing with sosfilt
def process_in_chunks(signal_data, chunk_size=100):
    # Design bandpass Butterworth filter (40–60 Hz)
    sos = signal.butter(4, [40, 60], btype='bandpass',
                        fs=sample_rate, output='sos')
    
    # Initialize filter state for sosfilt
    zi = signal.sosfilt_zi(sos)
    output = []
    
    for i in range(0, len(signal_data), chunk_size):
        chunk = signal_data[i:i+chunk_size]
        
        # If last chunk is shorter than chunk_size
        if len(chunk) < chunk_size:
            zi = zi[:, :chunk.shape[-1]]
        
        filtered_chunk, zi = signal.sosfilt(sos, chunk, zi=zi)
        output.extend(filtered_chunk)
    
    return np.array(output)

# Apply real-time simulation
real_time_filtered = process_in_chunks(signal_with_noise)

# Plot results
plt.figure(figsize=(12, 6))
plt.plot(t, signal_with_noise, label='Noisy Signal', alpha=0.5)
plt.plot(t, real_time_filtered, label='Real-Time Filtered (sosfilt)', linewidth=2)
plt.plot(t, clean_signal, '--', label='Original Clean Signal', alpha=0.7)
plt.xlabel("Time [s]")
plt.ylabel("Amplitude")
plt.title("Real-Time Bandpass Filtering Using sosfilt")
plt.legend()
plt.grid(True)
plt.tight_layout()
plt.show()

You can refer to the screenshot below to see the output:

iir filter python

This approach maintains the filter’s state between chunks, making it suitable for continuous processing of streaming data.

Read Python Scipy Odeint

Method 4 – Filter Analysis and Visualization

Understanding your filter’s characteristics is crucial. Let’s visualize the frequency response:

# Design a filter
sos = signal.butter(4, [40, 60], btype='bandpass', 
                   fs=sample_rate, output='sos')

# Calculate frequency response
w, h = signal.sosfreqz(sos, worN=2000, fs=sample_rate)

# Plot the frequency response
plt.figure(figsize=(10, 6))
plt.plot(w, 20 * np.log10(abs(h)))
plt.title('Filter Frequency Response')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Gain (dB)')
plt.grid(True)
plt.axvline(40, color='r', linestyle='--', alpha=0.7)
plt.axvline(60, color='r', linestyle='--', alpha=0.7)
plt.ylim(-60, 5)

This visualization helps you verify if your filter is behaving as expected in the frequency domain.

Read Python SciPy Interpolate

Method 5 – Cascade Multiple Filters

Sometimes a single filter can’t achieve the desired frequency response. In such cases, cascading multiple filters can be an effective solution:

# Design a highpass filter to remove baseline wander
sos1 = signal.butter(2, 0.5, btype='highpass', fs=500, output='sos')

# Design a lowpass filter to remove high-frequency noise
sos2 = signal.butter(2, 40, btype='lowpass', fs=500, output='sos')

# Design a notch filter to remove 60 Hz power line interference
sos3 = signal.iirnotch(60, 30, fs=500)

# Apply filters in sequence
filtered1 = signal.sosfilt(sos1, ecg_noisy)
filtered2 = signal.sosfilt(sos2, filtered1)
filtered3 = signal.sosfilt(sos3, filtered2)

# Plot
plt.figure(figsize=(10, 6))
plt.plot(t, ecg_noisy, 'b-', alpha=0.5, label='Noisy ECG')
plt.plot(t, filtered3, 'r-', label='Cascaded Filtering')
plt.legend()
plt.grid(True)
plt.title('ECG Filtering with Cascaded Filters')

I find this approach particularly useful when dealing with biomedical signals that have multiple noise sources at different frequencies.

Method 6 – Custom IIR Filter Design

While SciPy’s built-in filter types are sufficient for most applications, sometimes you may need to design a custom filter:

# Define desired frequency points
freqs = [0, 10, 20, 50, 100, 200, 500]  # Hz
gains = [-60, -30, 5, 5, -10, -60, -80]  # dB

# Convert to linear scale
gains_linear = 10**(np.array(gains)/20)

# Normalize frequencies to Nyquist
nyq = 500 / 2
freqs_norm = np.array(freqs) / nyq

# Design custom filter
b, a = signal.iirdesign(freqs_norm[2:4], 
                       [freqs_norm[1], freqs_norm[4]], 
                       gpass=3, gstop=40, 
                       ftype='cheby2')

# Convert to SOS format for better numerical stability
sos = signal.tf2sos(b, a)

# Apply filter
custom_filtered = signal.sosfilt(sos, signal_with_noise)

This example creates a custom filter with specific gain requirements at different frequency points.

Read Python SciPy Optimize Root

Method 7 – Forward-Backward Filtering for Zero-Phase

For offline processing where causality isn’t required, forward-backward filtering removes phase distortion:

# Design filter
b, a = signal.butter(4, [40, 60], btype='bandpass', fs=sample_rate)

# Regular filtering (with phase distortion)
forward_filtered = signal.lfilter(b, a, signal_with_noise)

# Forward-backward filtering (zero phase)
fb_filtered = signal.filtfilt(b, a, signal_with_noise)

# Plot to compare
plt.figure(figsize=(10, 6))
plt.plot(t, signal_with_noise, 'b-', alpha=0.3, label='Noisy Signal')
plt.plot(t, forward_filtered, 'g-', alpha=0.7, label='Forward Filtering')
plt.plot(t, fb_filtered, 'r-', label='Forward-Backward Filtering')
plt.plot(t, clean_signal, 'k--', alpha=0.5, label='Original Signal')
plt.legend()
plt.title('Phase Distortion Comparison')
plt.grid(True)

The difference is especially noticeable at signal transitions, where the forward-only filtered signal shows a delay.

Check out Python SciPy Fcluster

Practical Application: Clean ECG Signal

Let’s apply our knowledge to a real-world example – cleaning an ECG (electrocardiogram) signal from noise:

# Generate synthetic ECG-like signal (simplified)
def generate_ecg(duration, sample_rate):
    t = np.linspace(0, duration, int(duration * sample_rate), endpoint=False)
    # Main ECG wave components
    ecg = np.zeros_like(t)

    # Add periodic ECG peaks (simulating heartbeats)
    heart_rate = 70  # beats per minute
    beat_interval = 60/heart_rate

    for beat_time in np.arange(0, duration, beat_interval):
        # Find the closest time index
        idx = int(beat_time * sample_rate)
        if idx + 100 < len(ecg):
            # Add QRS complex
            ecg[idx:idx+10] += 0.1 * np.sin(np.linspace(0, np.pi, 10))  # Q wave
            ecg[idx+10:idx+30] += 1.0 * np.sin(np.linspace(0, np.pi, 20))  # R wave
            ecg[idx+30:idx+50] += 0.3 * np.sin(np.linspace(0, np.pi, 20))  # S wave
            ecg[idx+50:idx+100] += 0.2 * np.sin(np.linspace(0, np.pi, 50))  # T wave

    return t, ecg

# Generate the ECG signal
t, ecg_clean = generate_ecg(duration=5, sample_rate=500)

# Add different types of noise
power_line = 0.2 * np.sin(2 * np.pi * 60 * t)  # 60 Hz power line interference
baseline_wander = 0.5 * np.sin(2 * np.pi * 0.3 * t)  # Breathing artifacts
white_noise = 0.1 * np.random.normal(0, 1, len(t))  # Random noise

# Create noisy ECG
ecg_noisy = ecg_clean + power_line + baseline_wander + white_noise

# Design a bandpass filter to keep ECG frequencies (0.5-40Hz)
sos = signal.butter(4, [0.5, 40], btype='bandpass', fs=500, output='sos')
ecg_filtered = signal.sosfilt(sos, ecg_noisy)

# Plot
plt.figure(figsize=(12, 8))

plt.subplot(3, 1, 1)
plt.plot(t, ecg_clean)
plt.title('Clean ECG Signal')
plt.grid(True)

plt.subplot(3, 1, 2)
plt.plot(t, ecg_noisy)
plt.title('Noisy ECG Signal')
plt.grid(True)

plt.subplot(3, 1, 3)
plt.plot(t, ecg_filtered)
plt.title('Filtered ECG Signal')
plt.grid(True)

plt.tight_layout()

This example demonstrates how IIR filters can be used to remove common noise sources from medical signals like ECG data.

Check out Python SciPy Spatial Distance Cdist

Compare IIR and FIR Filters

It’s worth understanding when to use IIR versus FIR filters:

# Design comparable IIR and FIR filters
# IIR filter (4th order Butterworth)
sos_iir = signal.butter(4, [40, 60], btype='bandpass', fs=sample_rate, output='sos')

# FIR filter (with equivalent performance)
numtaps = 101  # FIR filters typically need many more coefficients
b_fir = signal.firwin(numtaps, [40, 60], pass_zero=False, fs=sample_rate)

# Apply filters
iir_filtered = signal.sosfilt(sos_iir, signal_with_noise)
fir_filtered = signal.lfilter(b_fir, 1, signal_with_noise)

# Compare computational efficiency
import time

# Time IIR filtering
start = time.time()
for _ in range(100):
    signal.sosfilt(sos_iir, signal_with_noise)
iir_time = time.time() - start

# Time FIR filtering
start = time.time()
for _ in range(100):
    signal.lfilter(b_fir, 1, signal_with_noise)
fir_time = time.time() - start

print(f"IIR filter time: {iir_time:.4f} seconds")
print(f"FIR filter time: {fir_time:.4f} seconds")
print(f"IIR is {fir_time/iir_time:.2f}x faster")

From my experience, IIR filters are significantly more efficient when steep cutoffs are needed, while FIR filters offer linear phase and guaranteed stability.

Read Python SciPy Pairwise Distance

Real-World Example: Audio Processing

Let’s apply IIR filtering to a common audio processing task – creating an equalizer:

# This would typically read from a file
sample_rate = 44100
duration = 3
t = np.linspace(0, duration, int(sample_rate * duration), endpoint=False)

# Create a synthetic audio signal with multiple frequencies
audio = (0.5 * np.sin(2 * np.pi * 100 * t) +  # Bass
         0.3 * np.sin(2 * np.pi * 1000 * t) +  # Mid
         0.2 * np.sin(2 * np.pi * 5000 * t))    # Treble

# Three-band equalizer
def apply_eq(audio, bass_gain, mid_gain, treble_gain):
    # Bass filter (lowpass at 300 Hz)
    sos_bass = signal.butter(2, 300, btype='lowpass', fs=sample_rate, output='sos')
    bass = signal.sosfilt(sos_bass, audio) * bass_gain
    
    # Mid filter (bandpass from 300 Hz to 3000 Hz)
    sos_mid = signal.butter(2, [300, 3000], btype='bandpass', fs=sample_rate, output='sos')
    mid = signal.sosfilt(sos_mid, audio) * mid_gain
    
    # Treble filter (highpass at 3000 Hz)
    sos_treble = signal.butter(2, 3000, btype='highpass', fs=sample_rate, output='sos')
    treble = signal.sosfilt(sos_treble, audio) * treble_gain
    
    # Combine all bands
    return bass + mid + treble

# Apply different EQ settings
audio_eq1 = apply_eq(audio, 1.5, 1.0, 0.7)  # Bass boost
audio_eq2 = apply_eq(audio, 0.7, 1.0, 1.5)  # Treble boost

# Visualize in frequency domain
from scipy.fft import rfft, rfftfreq

def plot_spectrum(signal_data, fs, label):
    N = len(signal_data)
    yf = rfft(signal_data)
    xf = rfftfreq(N, 1/fs)
    plt.semilogx(xf, 2.0/N * np.abs(yf), label=label)

plt.figure(figsize=(10, 6))
plot_spectrum(audio, sample_rate, 'Original')
plot_spectrum(audio_eq1, sample_rate, 'Bass Boost')
plot_spectrum(audio_eq2, sample_rate, 'Treble Boost')
plt.grid(True)
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude')
plt.title('Equalizer Frequency Response')
plt.legend()
plt.xlim(20, 20000)  # Audio range

This example demonstrates how IIR filters can be used to create a basic audio equalizer, a common application in music production and audio processing.

Conclusion

SciPy’s IIR filters provide a powerful and flexible toolset for digital signal processing in Python. Whether you’re cleaning biomedical signals, processing audio, or analyzing sensor data, these filters can effectively remove unwanted frequencies while preserving your signal of interest.

I’ve shown you seven different methods for implementing and analyzing IIR filters, from basic implementations to advanced techniques like cascading and custom filter design. Each approach has its strengths, and the best choice depends on your specific application.

Remember that IIR filters offer computational efficiency but come with trade-offs in phase response and potential instability. When in doubt, start simple and refine your approach based on the results.

You may also read:

51 Python Programs

51 PYTHON PROGRAMS PDF FREE

Download a FREE PDF (112 Pages) Containing 51 Useful Python Programs.

pyython developer roadmap

Aspiring to be a Python developer?

Download a FREE PDF on how to become a Python developer.

Let’s be friends

Be the first to know about sales and special discounts.