「包括的な時系列探索的分析」
「美容とファッションのトレンドを包括的に探る時系列分析」
大気品質データの詳細な分析
時間スタンプでインデックス化されたデータをお持ちです。データはストレージの需要と供給に関するものかもしれません。戦略的な製品の理想的な再補充間隔を予測するように指示されているかもしれません。または、歴史的な販売情報をチームに対してキープアクション可能な洞察に変換する必要があるかもしれません。おそらくデータは金融情報で、歴史的な金利情報と一部の株価の選択に関する情報を持っているかもしれません。市場のボラティリティをモデリングするよう指示されていて、投資期間にわたる通貨リスクを定量化する必要があるかもしれません。社会科学やエネルギー分布から医療まで、環境研究まで。例はたくさんあります。しかし、これらのシナリオには何が共通しているのでしょうか？第一に、時間の系列課題があります。そして、おそらく簡潔かつ包括的な探索的分析から始めることで確実に利益を得ることができます。
目次
この記事の目標
しかし、探索的な時系列分析を実行するとはどういうことでしょうか？他のデータサイエンスの問題とは異なり、時系列データからの洞察を得ることは難しく、簡単ではありません。データには重要な潜在的なトレンドや季節があるかもしれず、複雑な周期パターン内でのネストされた予測に適しているかもしれません。データ生成プロセスの障害による異常な外れ値と、重要な情報を持つ実際の異常値を区別することは難しいでしょう。また、欠損値の扱いも思ったほど簡単ではないかもしれません。
本記事では、私が時系列データセットを研究する際にうまく機能したプロセスを紹介します。私が詳細な視覚化と統計的要約を生成するための明確で高度な情報を提供するいくつかのベストプラクティスに焦点を当てます。
データセットの説明
ここで研究されたデータは、カナダのブリティッシュコロンビア州バンクーバー市の4つのモニタリングステーションからのものです。これらのデータには、2016年1月1日から2022年7月3日までの時間ごとの平均粒子状物質 PM 2.5（直径2.5マイクロメートル以下の微粒子）のµg/m3（立方メートル当たりのマイクログラム）の値が含まれています。
- 「Pythonを用いた巡回セールスマン問題の実装、解決、および可視化」
- このAI研究では、優れた画像分類のためにランダムスライスデータ拡張（RSMDA）を提案します：ニューラルネットワークの精度と堅牢性を向上させるための新しいアプローチ
- 「CNNによる特徴抽出の探求」
PM 2.5の主要な原因は化石燃料の燃焼であり、都市では通常、車の交通と建設現場から発生します。汚染物質の別の主要な源は森林火災や草原火災であり、風によって簡単に運ばれます [1]。
下の画像は、私たちが探索する駅のおおよその位置を示しています。
このデータセットはブリティッシュコロンビアデータカタログから取得されており、出版者によれば品質保証はされていません [5]。ここで表示されるバージョンでは、ネガティブな測定値（57kの観測値のうち6つのみ）を欠損値として割り当て、選択した駅のマスターデータフレームを生成するなど、いくつかの問題を前処理しています。
ライブラリと依存関係
視覚化のためにPython 3.9とプロットライブラリMatplotlibとSeabornを使用します。統計テストとデータ探索にはstatsmodels PythonモジュールとSciPyライブラリを使用します。データ操作と補助タスクにはPandasとNumpyを使用します。
これらのパッケージは、AnacondaやMiniconda、Google Collab、Kaggle Notebooksなどの人気のあるPythonディストリビューションとホストされたノートブックでネイティブに利用できます。したがって、ここにあるすべてのコード例は、選択した環境で簡単に再現できるはずです。
はじめに
ライブラリをインポートすることから始めましょう。matplotlib.dates.mdatesとdatetimeモジュールを呼び出し、DateTimeインデックスを操作するために使用します。一貫した視覚化を生成するために、プロットスタイルとカラーパレットを定義することもおすすめです。では、始めましょう。
.csvファイルを読み込んだ後、タイムスタンプのDATE_PSをNumPyのdatetime64オブジェクトとして定義し、DataFrameのインデックスに設定します。この一般的な手順により、Pandasの時系列機能が利用できるようになります。これを使用して、データセットの日付パートの特徴を作成するために以下で使用される手順などが可能になります。
全体像
ここで、私たちが深入りする場所の広範なビューがここにあります。これはデータを最初に把握するのに時間を費やすところです。
この可視化には、Seabornの関係プロットが使用されます。関係プロットは、DataFrameの集計された長いバージョンから読み取ります。そのために、Pandasのmeltとresampleメソッドを使用し、24時間ごとの平均集計でデータの粒度を時間ごとから日次の平均測定値に減らします。これはプロットの生成にかかる時間を短縮するために行われます。
全ての4つのステーションの時間範囲にわたる明確な図を持つことで、既にいくつかの観察を始めることが可能です：
- いくつかの重大な異常があり、それらは夏と初秋に普及しているようです。
- これらの異常は、大規模なイベントによるもののようで、それらはほぼ同じ期間に4つのステーションに影響を与えているようです。
- 注意して見てみると、各チャート内のすべてのステーションのグレーアウトした散布図が含まれています。この微細な詳細により、例えば2017年に存在する異常は、North Vancouverの両方のステーションに大きな影響を与えました（より高いPM 2.5値に達しました）、一方2018年のイベントでは逆のことが起こりました。このテクニックはまた、4つの折れ線グラフが同じY軸範囲内にあることも保証しています。
この最初のプロットから得られるいくつかの良い方法：
- Matplotlibは、高度にカスタマイズ可能な軸ティッカーロケーターとフォーマッターを許可します。この例では、
mdatesモジュールからMonthLocatorを使用して、X軸にマイナーマンスロケーターを作成し、全体の可読性を向上させています。 - プロットされた期間の正確な予定日範囲は、プロットのタイトルやサブタイトルに渡されます。これは可視化が切り詰められた場合に、プロットされた期間の両端のデータの欠落によるものかを確認するのに役立ちます。また、あなたの発見を報告するための良いドキュメンテーションの方法でもあります。
素晴らしいスタートですが、視点を少し絞りましょう。次のセクションでは、時間の短い期間を見ていきますが、元の時間単位（時単位）の粒度を持っています。
詳細な表示
これから、テーラーメイドの可視化を生成するために素早く呼び出すことができるいくつかの関数を定義し始めます。これは、今後非常に役立つ分析ツールセットをセットアップする手段と考えることができます。
まず最初の関数は、特定の時間範囲内の個々の時系列を見るのに役立ちます。North Vancouver Mahon Parkステーションの2017年を見てみましょう。
ここで、主要な異常に加えて、局所的な小さなスパイクも見られます。また、2017年の初めと12月には、より高い変動性（短い期間で分散が増加すること）の期間もありました。
異常範囲の外側を探求し、Y軸の狭い範囲内で値を分析できるようにしましょう。
ここではいくつかの欠落値が見られます。私たちの関数のfill=Trueパラメーターは、それを識別するのに役立ち、データの欠落性に視覚的強調を与える良い方法です。最初は気づきにくい小さな中断がはっきりと見えるようになりました。
また、X軸のカスタムな日付形式にも注目しているかもしれません。上記のプロットでは、カスタムなメジャーやマイナーロケーターおよびフォーマッターを使用して、私たちのplot_sequence()関数を高度化しました。この新機能により、プロットは可視化の期間に合わせて適応し、X軸が適切にフォーマットされます。以下に、関数に含まれていたコード断片を示します。
今、データセットには中断があることがわかっているので、それを詳しく見てみましょう。
欠損値
表形式のデータ問題では、このセクションではおそらくMAR（Missing At Random）とMNAR（Missing Not At Random）を定義することに焦点を当てることになるでしょう。ただし、データの性質（感覚的な時間的測定）を知っているので、データストリームの中断はおそらく意図的ではないとわかっています。したがって、こうした場合では、単独の欠損値と連続的な欠損値、完全に欠損しているサンプルとの区別が重要です。ここでは可能性が広がっていますので、将来的にはそれについての記事を専門に取り上げる予定です。
さて、まずは欠損値のヒートマップを見てみましょう。再度、そのための関数を定義します。
ヒートマップは欠損の程度を定量化し、時系列の中でそれらを特定するのに非常に役立ちます。ここから、以下のように記述できます:
- 完全に欠損しているサンプルはありません（欠損値が同時に一定期間に発生することはありません）。これは監視ステーションからのデータが独立して発生するため予想されます。
- 時間軸の早い段階では長期間の欠損値の連続があり、時間が経つにつれてデータの利用可能性が改善されているようです。
後のセクションでの一部の統計分析では、欠損値が問題となることがあります。そのため、Pandasのffill()やbfill()メソッドなど、必要に応じて単純な手法を使用してそれらを適切に処理することにします。
- 最も近い利用可能な値をそれぞれ前方または後方に伝達するために使用されるPandasの
ffill()およびbfill()メソッド。 - Pandasの
interpolate()メソッドを使用した線形またはスプライン補間。欠損した間隔を埋めるために隣接する観測値を使用して曲線を描画します。
断続性
データの性質からは、負の値は予想されないはずです。先ほど述べたように、データの前処理時にそれらを欠損として処理しました。さて、サマリー統計を呼び出して確認してみましょう。
各ステーションの最小測定値がゼロであることがわかります。これから次の問いに進みます。私たちの時系列は断続的なのでしょうか？
断続性は、データにゼロと完全に一致する値が多数含まれる場合に特徴付けられます。この振る舞いは特定の課題を引き起こし、モデル選択時に考慮する必要があります。では、私たちのデータにおいてゼロの値がどれだけ見られるのでしょうか？
ゼロの数がほとんどないことがわかりますので、私たちの時系列は断続的ではありません。これは容易ですが重要なチェックです。特に予測を行いたい場合には、予測モデルが絶対ゼロを予測することは困難かもしれません。たとえば需要を予測する場合、顧客に対して3つの製品を配送する計画を立てる必要がある場合に、実際には顧客が何も期待していない場合には問題が発生するかもしれません。
季節性
時系列のサイクルを理解することは、モデリングの段階を計画する上で基本的な要素です。データを集約しすぎると、より小さなサイクルの重要な情報を見失う可能性がありますし、より細かい粒度での予測の実現可能性を判断することもできます。
ここでは、いくつかのボックスプロットを使用して調査を開始します。ただし、まずは上位5%パーセンタイルを一時的に除外してデータをより良い尺度で見ることにします。
この次の関数では、データ内のサイクルを調査するために一連の箱ひげ図を使用します。また、中央値に対してカラーパレットをマッピングして、別の視覚的な手がかりとして使用します。
この最初のプロットでは、毎時の測定値が返されます。ここでは次のことが確認できます。
- 午前9時から午後2時までの間、PM 2.5の値が一貫して高くなっています。
- 北バンクーバー外のステーションでも、午後8時から午後10時にかけてピークがあります。
- 早朝には、午前2時から午前5時までの間にPM 2.5の最低値があります。
次に、週ごとの季節変動と週を通じた値の違いを見てみましょう。
ここから、次のことがわかります。
- 週末にはPM 2.5の値が低くなります。
- 火曜日には、汚染レベルがやや高くなる傾向があります。
最後に、月ごとのトレンドを見てみましょう。
ここでは次のことが観察されます。
- 8月には、すべての年で一貫して高いPM 2.5の値があります。
- 南部のステーションでは、6月と7月にPM 2.5の値が低くなりますが、北バンクーバーのステーションでは1月に測定値が低くなります。
最後に、これらのプロットからのより良い手法：
- カラーパレットを単純に
pallette="mako"として渡すだけでなく、その解釈を誤ってしまう可能性があるため、慎重に使用してください。カラーパレットを箱ひげ図に適用する際、X軸にマッピングされるのではなく、対象変数にマッピングされることに注意してください。 - グリッドプロットは、低次元データの情報を強力に保持する容器であり、Seabornの
relplot()やMatplotlibのsubplots()を使用して簡単に設定できます。 - Seabornの
boxplot()のorderパラメータを使用して、X軸を適切に並べ替えることもできます。私はそれを使って曜日のX軸ラベルを意味のある順序で再配置しました。
時間系列からトレンドの季節性をより詳細に把握するには、トレンド-季節分解を使用することができます。しかし、それは将来の記事で取り上げる予定であり、時間系列の類似性とモデル選択についてさらに深く探求できる場合に譲ります。興味があるトピックであれば、私の今後の記事を受け取るためにVoAGIで私をフォローしてください。
今のところ、4つのステーション間の線形関係を調査するために、よく知られた統計的係数のクイックビューを試してみましょう。
ピアソンの相関係数
Rのプログラマーは、次のプロットになじみがあるかもしれません。相関図は、GGallyパッケージのggpairs()のような複数のRライブラリで実装されている、簡潔で非常に情報量の多い視覚化手法です。相関図の上部の対角線には変数間の二変量相関、またはデータ内の数値変数間のピアソン相関係数が表示されます。下部の対角線には、データに適合した回帰曲線を持つ散布図が表示されます。最後に、主対角線には各変数のヒストグラムと密度曲線が表示されます。
以下のコードは、SeabornのPairGrid()プロットを使用した適応実装と、分析ツールセットの別の関数です。
予想通り、駅は強く相関しており、特にNorth Vancouverの2つの駅は近いため、相関が高いです。
計算時間を短縮するために、データは6時間ごとに集計されたことに注意してください。集計期間が大きいデータでこのプロットを試すと、データ中の外れ値を平滑化するために、相関係数が増加することがわかります。
時系列分析にすでに詳しい方は、他の種類の相関もチェックする価値があると思うかもしれません。ただし、まずは時系列の定常性をテストする必要があります。
定常性
定常系列とは、統計的な特性が時間とともに変わらない時系列のことを指します。つまり、定数の平均、分散、自己相関を持っているということです [4]。
多くの予測モデルは時系列の定常性に依存しているため、この探索的なフェーズでそれをテストすることが重要です。次の関数では、statsmodelsの実装を使用して定常性の2つの一般的に使用されるテスト、Augmented Dickey-Fuller（”ADF”）テストとKwiatkowski-Phillips-Schmidt-Shin（”KPSS”）テストを使用します。
以下に両テストの帰無仮説を示します。これらは対立する帰無仮説を持っているため、結果の容易な解釈のために「Decision」列を作成します。両テストについての詳細は、statsmodelsのドキュメントで読むことができます。
Augmented Dickey-Fuller（ADF）テストの帰無仮説:
• H0: 時系列サンプルに単位根が存在する（非定常）
• Ha: 時系列サンプルに単位根が存在しない（定常）
Kwiatkowski-Phillips-Schmidt-Shin（KPSS）テストの帰無仮説:
• H0: データは一定の周りで定常である（定常）
• Ha: 時系列サンプルに単位根が存在する（非定常）
今度は、どのスケールで定常性をチェックするかという重要な質問が浮かびます。その答えは、データのモデリング方法に大きく依存します。包括的な探索的分析の目標の1つは、この意思決定を支援することです。
説明のため、次の例では、バンクーバー国際空港駅の2016年1月と2022年1月のデータを見て、2016年から2022年にかけての振る舞いの変化があったかどうかを確認します。
Missing Valuesのセクションで説明したように、Pandasのffill()、bfill()、interpolate()メソッドを使用して、シリーズ内の欠損を迅速に補完することができます。両テストは完全なサンプルのみ受け入れるため、私は関数に専用の引数fillnaを定義して、これらのメソッドのどちらかを選択できるようにしました。
今、結果に戻ります。
2016年には両方のテストが非定常性を示していたが、2022年には結果が分かれた。 statsmodelsのドキュメントには、ADFテストとKPSSテストを同時に実行した場合の結果の解釈が明確に記載されている[6]：
• ケース1: 両方のテストが系列が非定常ではないと結論づける — 系列は非定常である
• ケース2: 両方のテストが系列が定常であると結論づける — 系列は定常である
• ケース3: KPSSは定常性を示し、ADFは非定常性を示す — 系列はトレンド定常である。 トレンドを除去することで系列を厳密に定常にする必要があります。トレンドの除去後の系列が定常性を持つかどうかを確認します。
• ケース4: KPSSは非定常性を示し、ADFは定常性を示す — 系列は差分定常です。系列を定常にするために差分を使用する必要があります。差分系列が定常性を持つかどうかを確認します。
この操作を複数か月にわたって4つの駅すべてに対して繰り返すと、ケース4がデータで優勢であることがわかります。これは、次のセクションで最初の階差を取ることでデータを定常化することにつながります。
最初の階差
最も一般的な変換技術の一つである、時系列に対して一次または二次の階差を適用することは、定常でない時系列に対してのみ使用できる統計モデルにデータを適用可能にするために広く使用されています。ここでは、2016年1月の前の例の1つに適用される技術を見てみましょう。ただし、最初に、変換前の元のデータをplot_sequence()関数で表示してみましょう。
月の始めから終わりまでの期間で分散が大きく変動していることがわかります。平均PM 2.5も高い値から低い、より安定した値に移動しているように見えます。これらは、系列の非定常性を確認するいくつかの特性です。
また、Pandasにはデータを差分化するための便利なメソッドがあります。DataFrameに.diff()を呼び出し、データの一次差分化バージョンとしてインスタンス化します。したがって、同じ期間を再度プロットしてみましょう。
まだ振動する分散ですが、データは明らかに平均値の周りで安定しています。差分化されたデータに対しても、stationarity_test()関数を再び呼び出して、定常性を確認することができます。
以上です。包括的な探索的時系列分析に対してさらなるチェックを加えました。次のことが確認できます：
- 非定常時系列を取り扱っています。
- 一次階差はそれを定常にする適切な変換技術です。
そして、最後のセクションに進むことにします。
自己相関
データが定常状態になったら、他の主要な時系列属性である「偏自己相関」と「自己相関」を調査することができます。と<。。 AC)は、時系列のラグ値間の線形関係を測定します。つまり、時系列自体との相関を測定します。[2]
偏自己相関関数（PACF）は、時間系列内のラグ値間の相関を測定します。ただし、それらの間で相関のある遅れた値の影響を除外します。これらは混乱変数として知られています。 [3]
両方の尺度は、相関グラフとして知られる統計プロットで視覚化することができます。ただし、まず、これらについての理解を深めることが重要です。
この記事は探索的分析に焦点を当てているため、これらの概念は統計予測モデルを構築する際に確かな直感を持つための非常に重要な考え方です。詳細な読解については、KaggleのノートブックグランドマスターLeonie Monigattiの素晴らしいカーネル「タイムシリーズ：ACFとPACFの解釈」をおすすめします。
上記の通り、自己相関は、時系列が以前のq個のラグとどのように相関するかを測定します。そのラグと自身のコピーがq期間だけ遅れるデータの部分集合との線形関係の測定と考えることができます。ACFは、移動平均（MA）モデルの次数qを決定するための重要な尺度です。 on>に対して、偏自己相関は時系列とそのpラグバージョンとの相関ですが、直接の影響だけに関係しています。たとえば、t-3からt-1の期間の偏自己相関を現在のt0の値とチェックしたい場合、t-3がt-2とt-1にどのように影響するか、t-2がt-1にどのように影響するかは気にしません。現在のタイムスタンプt0におけるt-3、t-2、およびt-1の直接の影響に焦点を当てます。部分自己相関（PACF）は、自己回帰（AR）モデルの次数pを決定するための重要な尺度です。
これらの概念が解説されたので、今度はデータに戻ることができます。これらの2つの尺度はよく一緒に分析されるため、最後の関数はPACFとACFのプロットを1つのグリッドプロットに統合し、複数の変数に対する相関グラフを返します。これには、statsmodelsのplot_pacf()とplot_acf()関数を使用し、それらをMatplotlibのsubplots()グリッドにマッピングします。
両方のstatsmodels関数は同じ引数を使用することに注意してくださいが、plot_pacf()プロット専用のmethodパラメーターは除外されています。
さて、データの異なる集約方法を実験してみることができますが、時間系列を再サンプリングすると、各ラグは異なる時間のジャンプを表すようになります。示しのために、2016年1月の4つのステーションのPACFとACFを、6時間の集約データセットで分析してみましょう。
相関グラフは、-1.0から1.0までの相関係数と、有意性のしきい値を示す影付きエリアを返します。それを超える値は統計的に有意であると見なすべきです。
上記の結果から、6時間間隔では次のことを結論づけることができます:
- ラグ1、2、3（t-6h、t-12h、およびt-18h）および時折4（t-24h）には有意なPACFがあります。
- ラグ1および4（t-6hおよびt-24h）はほとんどの場合において有意なACFを示しています。
そして、いくつかの最終的な良い慣行に注意してください：
- 高い粒度で長期時系列の相関図をプロットすることは避けるべきです（たとえば、時間ごとの測定値を持つデータセットの全年の相関図をプロットする場合）。サンプルサイズが増えるにつれて有意性の閾値はゼロに狭まります。
- 各ラグが表す時間期間でX軸に注釈を付けるために、私たちの関数に
x_labelパラメータを定義しました。その情報がない相関図を見ることも一般的ですが、それに簡単にアクセスできることで結果の誤解を避けることができます。 - Statsmodelsの
plot_acf()とplot_pacf()のデフォルト値はプロットに0ラグの相関係数を含めるように設定されています。任意の数の自己相関は常に1ですので、私たちのプロットはzero=Falseパラメータを使って最初のラグから始めるように設定しました。これにより、実際に分析する必要のあるラグがわかりやすくなり、Y軸のスケールも改善されます。
これにより、私たちは時間系列を徹底的に調査しました。視覚化と分析関数のツールセットを持っているため、データの包括的な理解を得ることができます。また、時間系列データセットを調査する際のベストプラクティスと、高品質のプロットで簡潔かつ洗練された形でそれらを提示する方法も学びました。
この記事を楽しんだ方は？
VoAGIで私をフォローすると、データサイエンス、機械学習、可視化、データ分析に関する記事をもっと読むことができます。
また、私をLinkedInで見つけることもできます。そこではこれらの内容の短いバージョンを共有しています。
参考文献
[1] “Department of Health — Fine Particles (PM 2.5) Questions and Answers.” Accessed October 14, 2022. https://www.health.ny.gov/environmental/indoors/air/pmq_a.htm.
[2] Peixeiro, Marco. “3. Going on a Random Walk.” Essay. In Time Series Forecasting in Python, 30–58. O’Reilly Media, 2022.
[3] Peixeiro, Marco. “5. Modeling an Autoregressive Process.” Essay. In Time Series Forecasting in Python, 81–100. O’Reilly Media, 2022.
[4] Peixeiro, Marco. “8. Accounting for Seasonality.” Essay. In Time Series Forecasting in Python, 156–79. O’Reilly Media, 2022.
[5] Services, Ministry of Citizens’. “The BC Data Catalogue.” Province of British Columbia. Province of British Columbia, February 2, 2022. https://www2.gov.bc.ca/gov/content/data/bc-data-catalogue.
[6] “Stationarity and Detrending (ADF/KPSS).” statsmodels. Accessed October 17, 2022. https://www.statsmodels.org/dev/examples/notebooks/generated/stationarity_detrending_adf_kpss.html.
We will continue to update VoAGI; if you have any questions or suggestions, please contact us!
Was this article helpful?
93 out of 132 found this helpful
Related articles
