Instantly illustrate the predominant period in time series data using spectrum analysis

Turn 2D data into a diagram instantly using python's matplotlib The second in the series, "Illustrate in an instant using python's matplotlib". Anyway, for those who want to easily ** find an outstanding cycle with one-dimensional time series data as shown below **.

image

Using python,

-[Spectrum analysis](https://ja.wikipedia.org/wiki/%E5%91%A8%E6%B3%A2%E6%95%B0%E3%82%B9%E3%83%9A%E3 Perform calculation of% 82% AF% E3% 83% 88% E3% 83% AB) --Illustration of calculation results

Do all at once.

1. Code

Set the following classes in advance.

spectra.py


# coding:utf-8

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


class Spectra(object):
    def __init__(self, t, f, time_unit):
        """
         - t :Time axis value
         - f :Data value
         - time_unit :Time axis unit
         - Po :Power spectral density
        """
        assert t.size == f.size  #Make sure that the length of the time axis and the length of the data are the same
        assert np.unique(np.diff(t)).size == 1  #Make sure all time intervals are constant
        self.time_unit = time_unit   #Unit of time
        T = (t[1] - t[0]) * t.size
        self.period = 1.0 / (np.arange(t.size / 2)[1:] / T)

        #Calculate power spectral density
        f = f - np.average(f)         #Zero the average.
        F = fftpack.fft(f)                          #Fast Fourier transform
        self.Po = np.abs(F[1:(t.size // 2)]) ** 2 / T

    def draw_with_time(self, fsizex=8, fsizey=6, print_flg=True, threshold=1.0):
        #Plot the power spectral density over time on the horizontal axis
        fig, ax = plt.subplots(figsize=(fsizex, fsizey))   #Specifying the size of the figure
        ax.set_yscale('log')
        ax.set_xscale('log')
        ax.set_xlabel(self.time_unit)
        ax.set_ylabel("Power Spectrum Density")
        ax.plot(self.period, self.Po)
        if print_flg:   #Draw a line where the power spectral density value is greater than the threshold and describe the period value.
            dominant_periods = self.period[self.Po > threshold]
            print(dominant_periods, self.time_unit +
                  ' components are dominant!')
            for dominant_period in dominant_periods:
                plt.axvline(x=dominant_period, linewidth=0.5, color='k')
                ax.text(dominant_period, threshold,
                        str(round(dominant_period, 3)))

        return plt

As you can see, FFT (Discrete Fourier Transform) in chronological order using scipy's fft function % 83% 95% E3% 83% BC% E3% 83% AA% E3% 82% A8% E5% A4% 89% E6% 8F% 9B) is performed to calculate the power spectral density. The period at which the power spectral density shows a larger value than the surroundings is the period that is predominant in the time series.

2. Verification

Create your own time series data and use it to verify the above code.

#coding:utf-8

import numpy as np
import matplotlib.pyplot as plt
from spectra import Spectra

#Creating time series data
#Monthly data for 30 years(That is, the number of data is 360)
#Big year(=12 months)10 years gradual in addition to the cycle(=120 months)A small amount of noise at each time
N = 360
t = np.arange(0, N)
td = t * np.pi / 6.0
f = np.sin(td) + 35.0 + 0.2 * np.sin(td * 0.1) + np.random.randn(N) * 0.1

#Drawing the original time series
plt.figure(figsize=(20, 6))
plt.plot(t, f)
plt.xlim(0, N)
plt.xlabel('Month')
plt.show()

#Drawing of outstanding cycles
spectra = Spectra(t, f, 'Month')
plt = spectra.draw_with_time()
plt.show()

3. Result

The original time series data diagram looks like the one at the top of this post. And the diagram of the outstanding cycle is as follows.

[ 120. 12.] Month components are dominant! image

As you can see, it ’s an outstanding cycle.

It extracted wonderfully. If the values at both ends of the original time series data are significantly different, [Window function](https://ja.wikipedia.org/wiki/%E7%AA%93%E9%96%A2%] Process the data into a form suitable for spectral analysis, such as by applying E6% 95% B0).


References

-[Practical climate data analysis using UNIX / Windows / Macintosh-Three points that writers of reports and papers on climate science, meteorology, ocean science, etc. should know](https://www.amazon. co.jp/UNIX-Windows-Macintosh%E3%82%92%E4%BD%BF%E3%81%A3%E3%81%9F%E5%AE%9F%E8%B7%B5-%E6%B0 % 97% E5% 80% 99% E3% 83% 87% E3% 83% BC% E3% 82% BF% E8% A7% A3% E6% 9E% 90-% E6% B0% 97% E5% 80% 99% E5% AD% A6-% E6% B0% 97% E8% B1% A1% E5% AD% A6-% E6% B5% B7% E6% B4% 8B% E5% AD% A6% E3% 81% AA% E3% 81% A9% E3% 81% AE% E5% A0% B1% E5% 91% 8A% E6% 9B% B8-% E8% AB% 96% E6% 96% 87% E3% 82% 92 % E6% 9B% B8% E3% 81% 8F% E4% BA% BA% E3% 81% 8C% E7% 9F% A5% E3% 81% A3% E3% 81% A6% E3% 81% 8A% E3 % 81% 8D% E3% 81% 9F% E3% 81% 843% E3% 81% A4% E3% 81% AE% E3% 83% 9D% E3% 82% A4% E3% 83% B3% E3% 83 % 88-% E6% 9D% BE% E5% B1% B1 / dp / 4772241221 / ref = sr_1_1? Ie = UTF8 & qid = 14482492333 & sr = 8-1 & keywords =% E6% B0% 97% E5% 80% 99% E3% 83 % 87% E3% 83% BC% E3% 82% BF% E8% A7% A3% E6% 9E% 90)

Recommended Posts

Instantly illustrate the predominant period in time series data using spectrum analysis
[Understand in the shortest time] Python basics for data analysis
Graph time series data in Python using pandas and matplotlib
Python: Time Series Analysis: Preprocessing Time Series Data
What you should not do in the process of time series data analysis (including reflection)
Time series analysis 3 Preprocessing of time series data
I tried logistic regression analysis for the first time using Titanic data
Power of forecasting methods in time series data analysis Semi-optimization (SARIMA) [Memo]
Get time series data from k-db.com in Python
Shortening the analysis time of Openpose using sound
How to read time series data in PyTorch
Predict from various data in Python using Facebook Prophet, a time series prediction tool
Big data analysis using the data flow control framework Luigi
Python: Time Series Analysis
Data analysis using xarray
Data analysis using Python 0
RNN_LSTM1 Time series analysis
Time series analysis 1 Basics
Analyzing the life of technology with Qiita article data ~ Survival time analysis using content logs ~
How to calculate the sum or average of time series csv data in an instant
Challenge to future sales forecast: ② Time series analysis using PyFlux
How to generate exponential pulse time series data in python
I tried to illustrate the time and time in C language
Time series analysis related memo
Time series analysis part 4 VAR
Time series analysis Part 3 Forecast
[Python] Plot time series data
Time series analysis Part 1 Autocorrelation
Data analysis using python pandas
[Unexpectedly known? ] Introducing a real day in the data analysis department
The first time a programming beginner tried simple data analysis by programming
<Pandas> How to handle time series data in a pivot table
[Statistics] [Time series analysis] Plot the ARMA model and grasp the tendency.
Try using FireBase Cloud Firestore in Python for the time being