医学と統計(49)

情報統計研究所へのアクセスはここからお気軽に

自己相関係数について。
前回は移動平均と時差相関を取り上げました。今回は、自己相関について考えてみたいと思います。時差相関では、末梢血中好塩基球数（Ba）と血清ＩｇＥ値の２群間のタイムラグでしたが、ここでの
自己相関は心電図波形から計測したＲ-Ｒ間隔時間（msec：以下、ＲＲと言う）を例題としています。通常、健常者であれば呼気と吸気時でＲＲが異なります。
図１は健常者の通常の呼吸時におけるＲＲの頻度を幹葉図で示しました。
　図１　通常呼吸時ＲＲ（SAS-JMP）

そして、
強制的にゆっくりとした深い腹式呼吸（約１０回／分）をおこなった時のRRの頻度を図２に示しました。
　図２　強制呼吸時RR（SAS-JMP）

通常呼吸時と強制呼吸時の統計量は次の通りでした。
                          通常呼吸時        強制呼吸時
平均値　 ＝  915.65 msec       914.67 msec
標準偏差＝    41.52 msec          73.75 msec
最頻値　 ＝   920     msec        830      msec

ここでの、
通常呼吸時と強制呼吸時におけるRRの自己相関係数は、図３と図４のようになりました。
これをコレログラムと言います。
　図３　通常呼吸時のコレログラム（SAS-JMP）

　図４　強制呼吸時のコレログラム（SAS-JMP）

強制呼吸時のコレログラムは規則的に循環変動していることが分かります。

自己相関係数の計算は前回の時差相関係数と同じで、RRを一つずつずらして計算したものですので、MSエクセル関数の「SUMSQ」と「SUMPRODUCT」を使って、前回の時差相関係数と同様の計算をすれば良いでしょう。

「R」なら、acf ( time series data ) で図５の結果を得ることができます。
しかし、
acf ( time series data ) の”time series data” は時系列データオブジェクトでなければなりません。もし、
RRのデータが例えば、次の様に並んでいるだけなら、
　dat <-  c ( 860,930,890,911,1000,891,・・・・ )

これを、
例えば、t_dat <-  ts ( dat, frequency=1 ) としてデータオブジェクトを作成しておき、
  acf ( t_dat )

とすれば、図５のような自己相関プロットが得られます。
　図５　「R」による自己相関プロット図

SAS-JMP なら、「分析→モデル化→時系列分析」で求めることが出来ます。

RRは副交感神経機能検査などで用いられており、ここでの自己相関分析も一つの検査方法として考えてみて下さい。

 

医学と統計(48)

情報統計研究所（統計分析のご相談ご依頼）のアクセスはここから。

前回の末梢血中好塩基球数（ Ba ) と血清IgE値（ IgE ) の実測データを用いて、実際に、時差相関係数を求めてみましょう。
Ba と IgE の測定単位が違いますので、次により、標準化（正規化変換）を行います。すなわち、
　Ba ( Xi ) の mean±SD=125.35±32.62、
　IgE ( Yi ) の mean±SD=583.48±183.42、
　
から、
　( Xi-125.35 ) ／ 32.62 、( Yi-583.48 ) ／ 183.42

により標準化（正規化変換）したデータを図１に示します。
図1　標準化した時系列データ

図１のデータについて、図２のようにデータを一つずつずらせて
　ｒ[k]＝Ｓxy[k]／SQRT（ Ｓｘｘ[0] ・Ｓyy[0] ）

を計算しますと、時差[k]（k=Lag Time）の相関係数が求まります。。
図２　時差相関係数の求め方

ＢａとＩｇＥの時差相関係数は次の通りです。
　ｒ[0]＝0.537、ｒ[1]＝0.549、ｒ[2]＝0.552、・・・、ｒ[10]＝-0.149

Ｂａ とＩ ｇＥ の時差はLag=2 で最大であり、Ｂａ の増加と ＩｇＥ の上昇には２日の時差があるようです。
さて、ＭＳエクセルの関数を使っての計算が面倒なら・・・、「Ｒ」を使用すれば簡単に求める事が出来ます。
「Ｒ」では ｃｃｆ（） が用意されており、
　ｃｃｆ （ X , Y )

で図３ の出力結果が得られます。
図３　「Ｒ」による時差相関プロット
 

時系列データでは、単相関で相関が見られなくても、時差を与える事で相関関係や循環性などの情報が得られる事もあります。

次回は自己相関について、心拍（R-R間隔）を取り上げたいと思います。