共计 4085 个字符，预计需要花费 11 分钟才能阅读完成。
周期是数据中呈现反复模式所需的工夫长度。更具体地说，它是模式的一个残缺周期的持续时间。在这篇文章中，将介绍计算工夫序列周期的三种不同办法。
咱们应用 City of Ottawa 数据集，次要关注的是每天的服务呼叫数量。所以不须要对病房名称进行初始数据处理。Ottawa 数据集在渥太华市提供的数据门户网站上收费提供。
让咱们加载 2019-2022 年的这些数据，并将它们连接起来失去一个 df。
fromgoogle.colabimportdrive
drive.mount('/content/gdrive')
importpandasaspd
importmatplotlib.pyplotasplt
importseabornassns
importnumpyasnp
file_path='/content/gdrive/My Drive/Colab Notebooks/Data/SR-2019.xlsx'
records2019=pd.read_excel(file_path)#,encoding='utf16')
file_path='/content/gdrive/My Drive/Colab Notebooks/Data/SR-2020.xlsx'
records2020=pd.read_excel(file_path)#,encoding='utf16')
file_path='/content/gdrive/My Drive/Colab Notebooks/Data/2021_Monthly_Service_Requests_EN.xlsx'
records2021=pd.read_excel(file_path)#,encoding='utf16')
file_path='/content/gdrive/My Drive/Colab Notebooks/Data/2022_Monthly_Service_Requests.csv'
records2022=pd.read_csv(file_path)
records=pd.concat([records2019,records2020,records2021,records2022],axis=0)让咱们依据服务调用日期聚合这些数据，并失去一个简略的图。
records["DATE_RAISED"]=pd.to_datetime(records.DATE_RAISED)
record_by_date=records.groupby("DATE_RAISED")["TYPE"].count().sort_index()
record_by_date.plot(figsize= (25, 10))
plt.ylabel('Number of requests')
plt.grid(visible=True,which='both')
plt.figure()
record_by_date.iloc[100:130].plot(figsize= (25, 10))
plt.ylabel('Number of requests')
plt.grid(visible=True,which='both')填充缺失
让咱们检查一下咱们的数据是否蕴含了所有的日期。
start_date=record_by_date.index.min()
end_date=record_by_date.index.max()
# create a complete date range for the period of interest
date_range=pd.date_range(start=start_date, end=end_date, freq='D')
# compare the date range to the index of the time series
missing_dates=date_range[~date_range.isin(record_by_date.index)]
iflen(missing_dates) >0:
print("Missing dates:", missing_dates)
else:
print("No missing dates")正如所预期的那样，数据短少一些日期的值。让咱们用相邻日期的平均值填充这些值。
# Reindex to fill missing dates
idx=pd.date_range(start=record_by_date.index.min(), end=record_by_date.index.max(), freq='D')
record_by_date=record_by_date.reindex(idx, fill_value=0)
# Add missing dates with average of surrounding values
fordateinmissing_dates:
prev_date=date-pd.DateOffset(days=1)
next_date=date+pd.DateOffset(days=1)
prev_val=record_by_date.loc[prev_date] ifprev_dateinrecord_by_date.indexelsenp.nan
next_val=record_by_date.loc[next_date] ifnext_dateinrecord_by_date.indexelsenp.nan
avg_val=np.nanmean([prev_val, next_val])
record_by_date.loc[date] =avg_val这就是咱们要做的所有预处理了，在所有这些步骤之后，咱们尝试检测这个工夫序列的周期。一般来说，基于假日模式和个别的人类习惯，咱们心愿在数据中看到七天的周期，咱们来看看是不是有这样的后果。
0、目测
最简略的办法就是目测。这是一种主观的办法，而不是一种正式的或统计的办法，所以我把它作为咱们列表中的原始办法。
如果咱们看一下这张图的放大部分，咱们能够看到 7 天的周期。最低值呈现在 5 月 14 日、21 日和 28 日。但最高点仿佛不遵循这个模式。但在更大的范畴内，咱们依然能够说这个数据集的周期是 7 天。
上面咱们来正式的进行剖析：
1、自相干剖析
咱们将绘制工夫序列的自相干值。查看 acf 图中各种滞后值的峰值。与第一个显著峰值对应的滞后能够给出周期的预计。
对于这种状况，咱们看看 50 个滞后值，并应用 statmodels 包中的办法绘制 acf。
fromstatsmodels.graphics.tsaplotsimportplot_acf
fig, ax=plt.subplots(figsize=(14,7))
plot_acf(record_by_date.values.squeeze(), lags=50,ax=ax,title='Autocorrelation', use_vlines=True);
lags=list(range(51))
ax.set_xticks(lags);
ax.set_xticklabels(lags);从上图能够看出，在 7、1、21 等处有峰值。这证实了咱们的工夫序列有 7 天的周期。
2、疾速傅里叶变换
对工夫序列进行傅里叶变换，寻找主频重量。主频率的倒数能够作为周期的估计值。
傅里叶变换是一种数学运算，它把一个简单的信号分解成一组更简略的正弦和余弦波。傅里叶变换广泛应用于信号处理、通信、图像处理以及其余许多迷信和工程畛域。它容许咱们在频域中剖析和操作信号，这通常是一种比在时域中更天然和直观的了解和解决信号的办法。
fromscipy.fftimportfft
# Calculate the Fourier transform
yf=np.fft.fft(record_by_date)
xf=np.linspace(0.0, 1.0/(2.0), len(record_by_date)//2)
# Find the dominant frequency
# We have to drop the first element of the fft as it corresponds to the
# DC component or the average value of the signal
idx=np.argmax(np.abs(yf[1:len(record_by_date)//2]))
freq=xf[idx]
period=(1/freq)
print(f"The period of the time series is {period}")输入为：The period of the time series is 7.030927835051545。这与咱们应用 acf 和目视查看发现的每周周期类似。
3、周期图
周期图 Periodogram 是一个信号或序列的功率谱密度 (PSD) 图。换句话说它是一个显示信号中每个频率蕴含多少总功率的图表。周期图是通过计算信号的傅里叶变换的幅值平方失去的，罕用于信号处理和频谱剖析。在某种意义上，只是后面给出的基于 fft 的办法的扩大。
fromscipy.signalimportperiodogram
freq, power=periodogram(record_by_date)
period=1/freq[np.argmax(power)]
print(f"The period of the time series is {period}")
plt.plot(freq, power)
plt.xlabel('Frequency (Hz)')
plt.ylabel('Power spectral density')
plt.show()周期图能够分明地看出，信号的最高功率在 0.14，对应于 7 天的周期。
总结
本文，咱们介绍了寻找工夫序列周期的三种不同办法，通过应用这三种办法，咱们可能辨认信号的周期性，并应用常识进行确认。
https://avoid.overfit.cn/post/2ae6a3c1b9824defbd013aecd0a70635
作者：Shashindra Silva
