如何对分段线性或非线性波形进行去趋势处理?

发布于 2025-01-13 16:26:41 字数 3219 浏览 0 评论 0原文

我试图消除波形中存在的趋势,如下所示:

在此处输入图像描述

为此,我使用 scipy.signal.detrend() 如下:

autocorr = scipy.signal.detrend(autocorr)

但是我没有看到任何明显的趋势趋平。我得到以下内容:

在此处输入图像描述

我的目标是从波形中完全消除趋势。我还需要对其进行概括,以便它可以消除任何类型的波形的趋势 - 无论是线性的、分段线性的、多项式的等等。

您能建议一种方法来做到这一点吗?

注意:为了复制上述波形,您只需运行我用来生成它的以下代码:

#Loading Libraries
import warnings
warnings.filterwarnings("ignore")

import json
import sys, os
import numpy as np
import pandas as pd
import glob
import pickle

from statsmodels.tsa.stattools import adfuller, acf, pacf
from scipy.signal import find_peaks, square
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
import matplotlib.pyplot as plt

#Generating a function with Dual Seasonality:
def white_noise(mu, sigma, num_pts):
    """ Function to generate Gaussian Normal Noise
    Args:
        sigma: std value
        num_pts: no of points
        mu: mean value

    Returns:
        generated Gaussian Normal Noise
    """
    
    noise = np.random.normal(mu, sigma, num_pts)
    return noise

def signal_line_plot(input_signal: pd.Series, title: str = "", y_label: str = "Signal"):
    """ Function to plot a time series signal
    Args:
        input_signal: time series signal that you want to plot
        title: title on plot
        y_label: label of the signal being plotted
        
    Returns:
        signal plot
    """

    plt.plot(input_signal)
    plt.title(title)
    plt.ylabel(y_label)
    plt.show()

# Square with two periodicities of daily and weekly. With @15min sampling frequency it means 4*24=96 samples and 4*24*7=672 

t_week = np.linspace(1,480, 480)
t_weekend=np.linspace(1,192,192)
T=96 #Time Period
x_weekday = 10*square(2*np.pi*t_week/T, duty=0.7)+10 + white_noise(0, 1,480)
x_weekend = 2*square(2*np.pi*t_weekend/T, duty=0.7)+2 + white_noise(0,1,192)
x_daily_weekly = np.concatenate((x_weekday, x_weekend)) 
x_daily_weekly_long = np.concatenate((x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly))
signal_line_plot(x_daily_weekly_long)
signal_line_plot(x_daily_weekly_long[0:1000])

#Finding Autocorrelation & Lags for the signal [WHICH  THE FINAL PARAMETERS WHICH ARE TO BE PLOTTED]:

#Determining Autocorrelation & Lag values
import scipy.signal as signal
autocorr = signal.correlate(x_daily_weekly_long, x_daily_weekly_long, mode="same")

#Normalize the autocorr values (such that the hightest peak value is at 1)
autocorr = (autocorr-min(autocorr))/(max(autocorr)-min(autocorr))

lags = signal.correlation_lags(len(x_daily_weekly_long), len(x_daily_weekly_long), mode = "same")

#Visualization
f = plt.figure()
f.set_figwidth(40)
f.set_figheight(10)
plt.plot(lags, autocorr)

#DETRENDING:
autocorr = scipy.signal.detrend(autocorr)

#Visualization
f = plt.figure()
f.set_figwidth(40)
f.set_figheight(10)
plt.plot(lags, autocorr)

I'm trying to remove the trend present in the waveform which looks like the following:

enter image description here

For doing so, I use scipy.signal.detrend() as follows:

autocorr = scipy.signal.detrend(autocorr)

But I don't see any significant flattening in trend. I get the following:

enter image description here

My objective is to have the trend completely eliminated from the waveform. And I need to also generalize it so that it can detrend any kind of waveform - be it linear, piece-wise linear, polynomial, etc.

Can you please suggest a way to do the same?

Note: In order to replicate the above waveform, you can simply run the following code that I used to generate it:

#Loading Libraries
import warnings
warnings.filterwarnings("ignore")

import json
import sys, os
import numpy as np
import pandas as pd
import glob
import pickle

from statsmodels.tsa.stattools import adfuller, acf, pacf
from scipy.signal import find_peaks, square
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
import matplotlib.pyplot as plt

#Generating a function with Dual Seasonality:
def white_noise(mu, sigma, num_pts):
    """ Function to generate Gaussian Normal Noise
    Args:
        sigma: std value
        num_pts: no of points
        mu: mean value

    Returns:
        generated Gaussian Normal Noise
    """
    
    noise = np.random.normal(mu, sigma, num_pts)
    return noise

def signal_line_plot(input_signal: pd.Series, title: str = "", y_label: str = "Signal"):
    """ Function to plot a time series signal
    Args:
        input_signal: time series signal that you want to plot
        title: title on plot
        y_label: label of the signal being plotted
        
    Returns:
        signal plot
    """

    plt.plot(input_signal)
    plt.title(title)
    plt.ylabel(y_label)
    plt.show()

# Square with two periodicities of daily and weekly. With @15min sampling frequency it means 4*24=96 samples and 4*24*7=672 

t_week = np.linspace(1,480, 480)
t_weekend=np.linspace(1,192,192)
T=96 #Time Period
x_weekday = 10*square(2*np.pi*t_week/T, duty=0.7)+10 + white_noise(0, 1,480)
x_weekend = 2*square(2*np.pi*t_weekend/T, duty=0.7)+2 + white_noise(0,1,192)
x_daily_weekly = np.concatenate((x_weekday, x_weekend)) 
x_daily_weekly_long = np.concatenate((x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly))
signal_line_plot(x_daily_weekly_long)
signal_line_plot(x_daily_weekly_long[0:1000])

#Finding Autocorrelation & Lags for the signal [WHICH  THE FINAL PARAMETERS WHICH ARE TO BE PLOTTED]:

#Determining Autocorrelation & Lag values
import scipy.signal as signal
autocorr = signal.correlate(x_daily_weekly_long, x_daily_weekly_long, mode="same")

#Normalize the autocorr values (such that the hightest peak value is at 1)
autocorr = (autocorr-min(autocorr))/(max(autocorr)-min(autocorr))

lags = signal.correlation_lags(len(x_daily_weekly_long), len(x_daily_weekly_long), mode = "same")

#Visualization
f = plt.figure()
f.set_figwidth(40)
f.set_figheight(10)
plt.plot(lags, autocorr)

#DETRENDING:
autocorr = scipy.signal.detrend(autocorr)

#Visualization
f = plt.figure()
f.set_figwidth(40)
f.set_figheight(10)
plt.plot(lags, autocorr)

如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论, 你可以免费 注册 一个本站的账号。

评论(1

祁梦 2025-01-20 16:26:42

由于它是自相关,因此它始终是偶数;因此,在 lag=0 处设置断点去趋势应该可以帮助你达到目标。

另一种消除趋势的方法是使用高通滤波器;你可以通过两种方式做到这一点。棘手的是决定截止频率应该是多少。

这是一种可能的方法:

#Loading Libraries
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

#Generating a function with Dual Seasonality:
def white_noise(mu, sigma, num_pts):
    """ Function to generate Gaussian Normal Noise
    Args:
        sigma: std value
        num_pts: no of points
        mu: mean value

    Returns:
        generated Gaussian Normal Noise
    """
    
    noise = np.random.normal(mu, sigma, num_pts)
    return noise

# High-pass filter via discrete Fourier transform
# Drop all components from 0th to dropcomponent-th
def dft_highpass(x, dropcomponent):
    fx = np.fft.rfft(x)
    fx[:dropcomponent] = 0
    return np.fft.irfft(fx)


# Square with two periodicities of daily and weekly. With @15min sampling frequency it means 4*24=96 samples and 4*24*7=672 

t_week = np.linspace(1,480, 480)
t_weekend=np.linspace(1,192,192)
T=96 #Time Period
x_weekday = 10*signal.square(2*np.pi*t_week/T, duty=0.7)+10 + white_noise(0, 1,480)
x_weekend = 2*signal.square(2*np.pi*t_weekend/T, duty=0.7)+2 + white_noise(0,1,192)
x_daily_weekly = np.concatenate((x_weekday, x_weekend)) 
x_daily_weekly_long = np.concatenate((x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly))

#Finding Autocorrelation & Lags for the signal [WHICH  THE FINAL PARAMETERS WHICH ARE TO BE PLOTTED]:

#Determining Autocorrelation & Lag values
autocorr = signal.correlate(x_daily_weekly_long, x_daily_weekly_long, mode="same")

#Normalize the autocorr values (such that the hightest peak value is at 1)
autocorr = (autocorr-min(autocorr))/(max(autocorr)-min(autocorr))

lags = signal.correlation_lags(len(x_daily_weekly_long), len(x_daily_weekly_long), mode = "same")


# detrend w/ breakpoints
dautocorr = signal.detrend(autocorr, bp=len(lags)//2)

# detrend w/ high-pass filter
# use `filtfilt` to get zero-phase
b, a = signal.butter(1, 1e-3, 'high')
fautocorr = signal.filtfilt(b, a, autocorr)

# detrend with DFT HPF
rautocorr = dft_highpass(autocorr, len(autocorr) // 1000)

#Visualization
fig, ax = plt.subplots(3)

for i in range(3):
    ax[i].plot(lags, autocorr, label='orig')

ax[0].plot(lags, dautocorr, label='detrend w/ bp')
ax[1].plot(lags, fautocorr, label='HPF')
ax[2].plot(lags, rautocorr, label='DFT')

for i in range(3):
    ax[i].legend()
    ax[i].set_ylabel('autocorr')

ax[-1].set_xlabel('lags')

给出

在此处输入图像描述

Since it's an auto-correlation, it will always be even; so detrending with a breakpoint at lag=0 should get you part of the way there.

An alternative way to detrend is to use a high-pass filter; you could do this in two ways. What will be tricky is deciding what the cut-off frequency should be.

Here's a possible way to do this:

#Loading Libraries
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

#Generating a function with Dual Seasonality:
def white_noise(mu, sigma, num_pts):
    """ Function to generate Gaussian Normal Noise
    Args:
        sigma: std value
        num_pts: no of points
        mu: mean value

    Returns:
        generated Gaussian Normal Noise
    """
    
    noise = np.random.normal(mu, sigma, num_pts)
    return noise

# High-pass filter via discrete Fourier transform
# Drop all components from 0th to dropcomponent-th
def dft_highpass(x, dropcomponent):
    fx = np.fft.rfft(x)
    fx[:dropcomponent] = 0
    return np.fft.irfft(fx)


# Square with two periodicities of daily and weekly. With @15min sampling frequency it means 4*24=96 samples and 4*24*7=672 

t_week = np.linspace(1,480, 480)
t_weekend=np.linspace(1,192,192)
T=96 #Time Period
x_weekday = 10*signal.square(2*np.pi*t_week/T, duty=0.7)+10 + white_noise(0, 1,480)
x_weekend = 2*signal.square(2*np.pi*t_weekend/T, duty=0.7)+2 + white_noise(0,1,192)
x_daily_weekly = np.concatenate((x_weekday, x_weekend)) 
x_daily_weekly_long = np.concatenate((x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly))

#Finding Autocorrelation & Lags for the signal [WHICH  THE FINAL PARAMETERS WHICH ARE TO BE PLOTTED]:

#Determining Autocorrelation & Lag values
autocorr = signal.correlate(x_daily_weekly_long, x_daily_weekly_long, mode="same")

#Normalize the autocorr values (such that the hightest peak value is at 1)
autocorr = (autocorr-min(autocorr))/(max(autocorr)-min(autocorr))

lags = signal.correlation_lags(len(x_daily_weekly_long), len(x_daily_weekly_long), mode = "same")


# detrend w/ breakpoints
dautocorr = signal.detrend(autocorr, bp=len(lags)//2)

# detrend w/ high-pass filter
# use `filtfilt` to get zero-phase
b, a = signal.butter(1, 1e-3, 'high')
fautocorr = signal.filtfilt(b, a, autocorr)

# detrend with DFT HPF
rautocorr = dft_highpass(autocorr, len(autocorr) // 1000)

#Visualization
fig, ax = plt.subplots(3)

for i in range(3):
    ax[i].plot(lags, autocorr, label='orig')

ax[0].plot(lags, dautocorr, label='detrend w/ bp')
ax[1].plot(lags, fautocorr, label='HPF')
ax[2].plot(lags, rautocorr, label='DFT')

for i in range(3):
    ax[i].legend()
    ax[i].set_ylabel('autocorr')

ax[-1].set_xlabel('lags')

giving

enter image description here

~没有更多了~
我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站,即表示您同意使用 Cookies 和您的相关数据。
原文