2.9 Filtering Data

Introduction¶

Filtering is an essential technique in time series analysis, especially in geosciences, where signals are often contaminated by noise or contain both short- and long-duration components. Filters help isolate meaningful information from raw data by attenuating unwanted frequencies or enhancing certain features. In this lecture, we will filter time series data, with a focus on a climate variable that has seasonality and a positive trend. Climate data often exhibits both short-term variations (such as daily or seasonal cycles) and long-term trends (such as oceanic warming). By applying filters, we can focus on specific components of a climate-related signal, whether we are interested in long-term climate trends or short-term weather patterns.

Why Filtering is Important¶

• Noise reduction: Geophysical data often contain noise from measurement instruments, environmental conditions, or unrelated signals. Filtering helps remove this noise and enhances the signal of interest.
• Feature isolation: By focusing on specific frequency bands, filtering allows us to isolate short-term phenomena (like storms) or long-term processes (like climate trends).
• Smoothing data: In geosciences, smoothing noisy time series data makes patterns more apparent and improves the clarity of visualizations.
• Detecting trends: Long-term filtering can reveal underlying trends in the data, which matters for studies of climate change, ocean circulation, and global warming.

Types of Filters¶

• Low-pass filter: Allows low-frequency components to pass through while attenuating high-frequency components. Useful for isolating long-term trends.
• High-pass filter: Allows high-frequency components to pass through while attenuating low-frequency components. Useful for focusing on short-term fluctuations.
• Band-pass filter: Allows a specific range of frequencies to pass through, blocking both higher and lower frequencies. Useful for analyzing phenomena within a particular frequency range.
• Smoothing filters: Such as moving averages or Gaussian filters, smooth the data to remove short-term fluctuations.

The data may superimpose multiple signals of various frequencies. To remove or extract specific signals that do not overlap in frequencies, we can filter the data.

The filter can be:

• high pass: reduce signals at frequencies lower than a corner frequency $f_c$, only let the signals above $f_c$. Often parameterized in functions as hp or highpass.
• low pass: reduce signals at frequencies greater than a cutoff frequency $f_c$, only let the signals below that $f_c$. Often parameterized in functions as lp or lowpass.
• band pass: reduce signals at frequencies lower than a low corner frequency $f_{c1}$ and at frequencies greater than a high corner frequency $f_{c2}>f_{c1}$. Often parameterized as bp or bandpass.

There exist different types of filters. The most common are butterworth and chebyshev, but there exist others.

Figure: Examples of filters illustrated here.

Example 1: Filtering a Synthetic Climate Time Series¶

We build a synthetic daily climate series with three known components: a linear warming trend, a seasonal cycle, and noise. Because we build it ourselves, we know the true components exactly, so we can check how well filtering recovers each one.

We make the noise colored (red) rather than white: real climate noise has more power at low frequencies than at high frequencies. We generate it by shaping the spectrum of white noise in the Fourier domain — scale the amplitude spectrum by $1/f$, keep random phases, and inverse transform. This choice makes the exercise honest: red noise has power at all frequencies, including the low frequencies where the trend and seasonal cycle live, so no filter can separate the components perfectly.

🖥️ Lecture slides — Session 08 (Fri Oct 16)

Recovering each component with filters¶

The three components occupy different frequency bands (frequencies here are in cycles per day):

• the trend lives at very low frequencies, near 0;
• the seasonal cycle is a narrow line at $f = 1/365 \approx 0.0027$ cycles per day;
• the noise spreads across all frequencies.

So we design three Butterworth filters:

• a low-pass with cutoff below the seasonal frequency to recover the trend;
• a band-pass bracketing $1/365$ cycles per day to isolate the seasonal cycle;
• a high-pass with cutoff above the seasonal frequency to isolate the noise.

We use butter(..., output='sos'), which returns the filter as second-order sections: a product of second-order polynomials that represents the same filter but is numerically more stable than the (b, a) polynomial form, especially at high filter orders. We apply the filter with sosfiltfilt, which runs the filter forward and backward so the result has zero phase shift (more on that below).

Why the recovery is imperfect: spectral leakage between components¶

The recovered trend is not a straight line: it wanders around the true trend. The recovered noise is missing part of the true noise. This is not a bug in the filters — it is a property of the data.

A filter separates signals by frequency band. It can only separate components cleanly if they occupy disjoint bands. Here the red noise has power at all frequencies, including below the low-pass cutoff. That low-frequency part of the noise passes through the low-pass filter along with the trend, and there is no way for the filter to tell them apart. The same leakage contaminates the band-pass output: the recovered “seasonal” series contains the noise power that falls inside the pass band, so its amplitude fluctuates from year to year even though the true seasonal cycle is perfectly regular. Meanwhile the high-pass output misses the low-frequency part of the noise — and red noise holds most of its power at low frequencies, so the high-pass filter recovers only a small fraction of the true noise and its RMS error is nearly as large as the noise itself.

Also note the filter edges: near the start and end of the series, the filter has incomplete data and the output is distorted. Edge effects are a standard filtering artifact; taper or trim the edges before interpreting them.

The lesson: filtering recovers a frequency band, not a physical component. The two coincide only when the components are spectrally separated. With colored noise, some leakage is unavoidable, and you should report it rather than ignore it.

Zero-phase versus causal filtering¶

scipy.signal gives two ways to apply an SOS filter:

• sosfilt runs the filter forward in time only. This is a causal filter: the output at time $t$ depends only on samples up to $t$. This is the only option in real-time systems, and it preserves the onset of a signal — nothing appears in the output before it appears in the input. The cost is a frequency-dependent phase delay: features in the output arrive late.
• sosfiltfilt runs the filter forward, then backward. The two passes cancel each other’s phase delay, so the output has zero phase: filtered features stay aligned in time with the raw data. The cost is that the filter is no longer causal — energy leaks backward in time, so a sharp onset acquires a small precursor. Use it for offline analysis when timing alignment matters, never for real-time processing, and be careful when picking arrival times near sharp onsets.

We demonstrate the difference on the seismic data below, where the sharp P-wave onset makes the phase delay easy to see.

Use Cases in Geoscience:¶

• Climate change studies: Filtering temperature data to remove noise and focus on long-term trends.
• El Niño and La Niña detection: Filtering helps identify periodic oscillations in sea surface temperature data.
• Weather forecasting: High-pass filters isolate short-term variations for analysis.

Example 2: Application to Seismology¶

We download the same data as in lecture 2.8: seismograms recorded in the Puget Sound (station UW.RATT) for the M8.2 Chignik, Alaska earthquake of July 29, 2021, plus a noise window before the event. We query the FDSN data center with the IRIS client; note that IRIS data services are now operated by EarthScope, but the client name IRIS still works.

The channel code HHZ denotes a high-gain broadband seismometer, sampled at 100 samples per second, vertical component. We keep the data in raw digitizer counts here (we do not remove the instrument response), so amplitude axes are labeled in counts.

We will use the scipy.signal module to filter the time series.

We use a butterworth filter of second order, band-passed between the frequencies of 1 Hz and 10 Hz. The sos output is a second-order sections representation: the filter is expressed as a product of second-order polynomials, which is numerically more stable than the single high-order polynomial form. We apply the filter with sosfiltfilt for a zero-phase result.

Now filter in a higher frequency band (10-40 Hz) and compare earthquake and noise signals.

The earthquake stands out from the noise in both frequency bands. In the 10-40 Hz band, the P wave is the dominant arrival: high frequencies attenuate with distance, so the later, slower phases are depleted in high-frequency energy.

Zero-phase (sosfiltfilt) versus causal (sosfilt) filtering¶

We now apply the same 1-10 Hz band-pass filter both ways and zoom in on the P-wave onset.

In the zoomed panel, the causal sosfilt output is shifted later in time relative to the zero-phase sosfiltfilt output: the filter’s phase delay moves the apparent onset. The zero-phase version stays aligned with the raw data, but it achieves this by filtering backward in time, which smears a small amount of energy before the onset.

Practical guidance:

• Use causal filtering (sosfilt) for real-time applications (earthquake early warning) and when the presence or absence of energy before an onset matters. If you pick arrival times on causally filtered data, correct for the filter delay or accept a bias.
• Use zero-phase filtering (sosfiltfilt) for offline analysis where timing alignment between raw and filtered data matters. Do not interpret small precursory wiggles near sharp onsets: they can be filter artifacts.

Example 3: Filtering imperfect records — a gap and a clock error¶

Everything above assumed a continuous, correctly timed record. Real archives do not oblige: telemetry drops out, disks fill, GPS clocks lose lock. The RATT record happens to be complete and well timed — so we break a copy of it on purpose. That way the unbroken record is the ground truth, and every repair gets graded.

A gap in the record¶

We knock out 20 seconds of the coda and fill it with zeros, which is how many archives (and a careless merge) deliver dropouts. Then we filter straight across, as if nothing happened.

The damage is out of all proportion to the gap. Zero-filling creates two step discontinuities whose height is the raw amplitude — dominated here by the long-period surface waves, tens of thousands of counts — while the true 1–10 Hz signal in this part of the coda is under a hundred counts. A step is broadband, so the filter responds with its own ringing at the corner frequencies: within two seconds of each edge the output is more than a hundred times the signal, and inside the gap the naive output shows plausible-looking oscillations that are pure filter artifact. As in lesson 2.6, the fabricated samples are the dangerous ones — they look like data.

No filter can resurrect the 20 seconds that were never recorded; repair means confining the damage. The fix is segment-wise filtering: filter each contiguous stretch of real data separately, and leave the gap as NaN. sosfiltfilt pads each segment internally (odd reflection about the endpoints), so no segment ever sees a step.

Graded against the complete record: segment-wise filtering cuts the near-edge error by two orders of magnitude, and beyond two seconds from the gap it matches the truth to better than a fraction of a percent. The residual contamination spans roughly two seconds — a few periods of the lowest passband frequency (1 Hz), which is the memory of the filter. The practical recipe: filter segments, keep gaps as missing, and flag a guard interval of a few filter periods at each segment edge in your metadata, so downstream users (or your own feature-extraction code in 2.11) know which samples to trust.

A clock error¶

The second pathology leaves the waveform perfect and corrupts only the timestamps. GPS clock drift, leap-second bugs, and digitizer restarts routinely produce sub-second timing errors — invisible to the eye on any plot, and fatal for phase picks, cross-correlation tomography, and any ML label derived from arrival times. We simulate an instrument whose clock is half a second late, and recover the offset the standard way: cross-correlate against a reference and read the lag of the peak.

The cross-correlation peak sits at exactly the injected offset: timing errors that no visual inspection would catch are recoverable to the sample, and to a fraction of a sample if you interpolate the correlation peak (fit a parabola through the three points around the maximum). On real data the reference is a co-located sensor, a neighboring station corrected for the predicted travel-time difference, or — for ambient-noise methods — the long-term average of the noise correlation itself, which is how modern networks monitor clock drift continuously. The daily-cadence version of this pathology (a GNSS series reporting the field at the wrong day) is one of the injectable defects in mlgeo_synth.degrade_series, and the repair logic is the same: detect by correlation against a reference, then shift.

Conclusion¶

Filtering time series data is a routine step in geoscientific analysis. It isolates frequency bands, removes noise, and reveals trends — but it separates bands, not physical components, and every filter choice (corner frequencies, order, causal vs zero-phase) leaves a signature in the output. Know that signature before you interpret the result. And before filtering at all, check the record itself: filter across a zero-filled gap and the artifact will outweigh the signal; trust a timestamp without a reference and a half-second error will ride silently into every product downstream.