地磁気絶対観測における基線値の異常値判定基準の定量化に向けて
地磁気観測所テクニカルレポート
第１巻第１号
１９
頁 平成１５年３月      
１

１

１
９


２００３


１




伊藤信和（鹿屋出張所）・藤井郁子（調査課）
２００３年１月３１日受付，２００３年３月３日改訂，２００３年３月７日受理

要旨

地磁気絶対観測において観測基線値に含まれる異常値を客観的に検出するため，ロバスト統計学を適用して，異常値の判定基準を定量化することを試みた．使用したデータは，鹿屋出張所で２０００年２月２５日〜１２月３１日に観測された ，，成分の観測基線値データで，・成分は５８回分４６４個，成分は７２回分５７６個である．各回の観測の標本数は８個程度と少数でありロバスト統計を適用するのは難しかったが，１回の観測ごとに観測基線値から中央値を除いた残差を解析期間中にわたって集めた標本集団は，少数の標本を除いて正規分布的であることがわかった．そこで，約１年分の残差を使って異常値検出を行った．
約１年分の残差の プロット（プロット）を作成し，標本分布が異常値のため正規分布的でなくなる部分を調査した．プロットの連続性が途切れるところをしきい値と設定したところ，，成分で中央値からの差の大きさが中央絶対偏差の７倍，成分で５倍以上を異常とみなすことになり，異常と判定された標本数は 成分が８個，成分が１３個，成分が１２個であった．異常値を除いて計算した残差の平均と標準偏差を使って，異常値と平均との差を標準偏差で規格化してみると，本報告で用いた異常値判定は標準偏差の３倍と同程度かやや厳しい基準に基づいていることがわかった．
ロバスト推定の考え方を導入することにより，客観的に異常値を検出でき，また１回ごとでは埋もれていた異常値の検出も可能になることが示唆された．本報告で扱った２０００年の観測基線値
データに関しては，中央絶対偏差の５〜７倍が異常値判定のしきい値の目安になると考えられるが，さらに精度の良い見積もりのためにはデータ数を増やしてしきい値を確定する必要がある．

１．はじめに地磁気観測所では，ある時刻の地磁気ベクトルを観測するために，ベクトルの変化量を連続的に測る変化観測と，基線値を得ることを目的とした絶対観測を行っている．絶対観測はおよそ週１度の頻度で行われ，プロトン磁力計及び角度測定器を用いて地磁気ベクトル３成分の値を測定し，変化観測で得られた同時刻の地磁気と比較して観測基線値を導出している．
１回の絶対観測では，３成分ごとに原則として８個の観測基線値を平均することによって，その回の基線値を求めている．しかし，８個の観測基線値の中に何らかの原因で異常な値が含まれていると，真の基線値を推定するのに障害となる．そのため異常値を検出し取り除く操作（オミット）や再観測が取り入れられているが，異常値の検出に定量的な基準がないため観測者によるオミットの判定基準のばらつきが指摘されてきた．
本報告では，異常値が含まれるデータのためのロバスト推定（，１９８１）の考え方を実際の観測基線値データに適用し，異常値のオミット基準を定量化することを試みた．
２．ロバスト推定
観測値の集合に異質な性質をもった異常値が含まれているとき，異常値の影響を除いた観測値分布の性質を調べるロバスト統計学（，１９８１）が有用となる．本報告で適用した手法を，異常値が含まれる観測値分布のパラメーター，異常値の検出法，異常値の影響の除去法にわけて説明する．
２．１異常値が含まれる観測値分布のパラメー
ター
ある物理量 を推定するため，回測定をくり返２

伊藤信和・藤井郁子

し観測値（標本）の集合｛｝，
＝１〜 を得たとす
る．観測エラーは平均０，標準偏差 の正規分布に従うランダムエラーだけとすると，物理量 と測定のばらつきの尺度 の最も良い推定量は，それぞれ，標本平均μと標本標準偏差σで与えられる．ただし，
１ １ ，σ＝（）．
μ
＝


−μ
Σ

−
１Σ

＝
１


＝
１
しかし，観測エラーに の分布に従わない異常値が含まれている場合，μとσは異常値も演算の中に取り入れて計算するので，異常値の振幅によってはや からかけ離れた値になる可能性がある．
異常値の影響を受けにくい統計的な推定量をロバスト推定量という．上述のように，μとσはロバストな推定量ではない．一方，中央値と中央絶対偏差（：）は異常値の影響を受けにくく，それぞれ と のロバストな推定量であることが知られている．ここで，中央値 とを大きい順に並べたときの ２
番は，個の標本

目の値であり，中央絶対偏差は
と  の残差の絶対値｜
｜を０を除いて大きさ順に並べたとき
− の中央値と定義される．
 と がロバスト推定量であることを，例を使って示す．下に示した標本集合は鹿屋出張所において２０００年４月１９日の絶対観測で求められた 成分の観測基線値（単位：）である．
｛３３．
３３，３４．
６１，３４．
５４，３４．
６２，３４．
４１，３４．
６８，３４．
７９，３４．
５９｝
この８個の標本の平均値μ，標準偏差σ，中央値
はそれぞれ，３４．
４５，０．
４６，，中央絶対偏差 ３４．
６０，０．
０７である．次に，１個目の値３３．
３３は本報告の解析の結果，異常であるとみなされた値であるので（第４章参照），この値を除いて７個の平均値μ，標準偏差σ，中央値 ，中央絶対偏差 を求めると，それぞれ，３４．
６１，０．
１２，３４．
６１，０．
０７となる．
このとき平均値の変化率は０．
４６％であるが，中央値の変化率は０．
０３％である．この例より，中央値は平均値ほど異常値の影響を受けないことがわかる．
平均値は０．
１６変化したが，異常値を除去した後の標準偏差は０．
１２であることから，この例における平均値の変化は標本のばらつき具合に対して有意な量である．標準偏差と中央絶対偏差も同様で，σとσには約４倍の差があり標本のばらつき具合がまったく異なる分布の印象を与えるが，と を比べると差がなく，異常値による影響が小さいことがわかる．

次にこの標本集合において中央値が真の基線値を代表しているかどうかを考える．異常値を除去した標本集合が本来の分布に従ったものとすると，中央値は異常値除去後の平均値に近い値になることから，中央値はこのような異常値を含む標本集合の真の基線値の推定量としての役割を果たすと考えられる．
２．２ 異常値検出法
μやσを用いて標本に含まれる異常値を検出する
が正規分布に従う場方法がある．ある標本集団


がμ±
３σから外れる確率は０．
３％しかないこ合，とから，｜

μ
｜
＞
３σ
−

（１）

を満たす

を異常値とする方法（３σエディット・ルール）や，そのしきい値が標本数によって変化するトンプソンの棄却検定などである．しかし，前節で示したようにμとσは異常値のため標本の分布を正しく示さないことがあるので，これらの方法はしばしば破たんすることが指摘されている（，２００２）．
μとσのかわりにロバストな  と を使え
ば，異常値が含まれる標本集団でも標本の分布を客観 的・定 量 的 に 見 積 も る こ と が で き る．３σ エ
ディット・ルールを構成するμとσをそれぞれ とに置き換たものは，ハンペルの判定法と呼ばれる（，２００２）．本報告も同様の手法で解析を行っているが，しきい値を３に固定せず，
｜＞｜

−

（２）

を満たす
を異常値とした．はしきい値を決定するパラメーターである．
２．３ 異常値の影響の除去法
２．２節の方法で標本の中に異常値が検出された
とき，異常値の影響を除いて標本集合の統計的な性質を求める方法には，大きくわけて２種類ある．一つは，異常値を完全に取り除くことで標本が の
本来の分布に従うようにする 
 法である．現行のオミットはこの方式である．二つ目は，異常値に重みをつけることで標本分布を本来の分布に近付ける 
 法である．
本報告では異常値の除去法について特に検討しなかった．
３．観測基線値の性質
図１に，本報告で扱う鹿屋出張所の２０００年２月２５

地磁気絶対観測における基線値の異常値判定基準の定量化に向けて

図１ ２０００年２月２５日〜１２月３１日に鹿屋出張所で観測された H（上），D（中），Z（下）成分の観測基線値，中央値及び MAD の時間変化．グラフ上段には観測基線値を黒丸，中央値を白丸で，また下段には MAD を棒グ
ラフで示す．

３

４

伊藤信和・藤井郁子

日〜１２月３１日の観測基線値，それらの中央値，およたデータがあることがわかる．
び を示す．ここでは異常値と思われるデータ
一回の絶対観測の目的は，個の観測基線値
，
が存在することがわかっているので，ロバストな推＝
１，から真の基線値 を推定することである
は ，あるいは を指す）．前述のように，現定量である中央値と を用いる．地磁気は大き（在の方式では基線値の推定量μは 個の観測基線値さと方向を持つベクトル量であり，当観測所では水の平均によって求められている．個のうちに他の平成分（），鉛直成分（），偏角（）の３成分で表している．観測基線値の単位は，・成分は ， 観測基線値と大きく異なるデータが得られた場合は異常値とみなし，異常値を除いた平均によって基線成分は分である．
値を得る．
絶対観測はおおむね１週間ごとに行われ，原則と異常値の検出基準の定量化のため，まず，観測基して８個の観測基線値が得られるまで測定をくり返
線値がどのような統計分布に従っているのかを調べすが，所々で再観測が行われているため，その頻度た．図２はある平均的な一日の観測基線値が正規分や一回の観測での観測基線値の個数は一定ではな）
布に従う仮定のもとに，確率密度関数 （い．図１の期間に絶対観測は ・成分で５８回，成
分で７２回行われ，得られた観測基線値の総数は，１（μ）
−）＝
（−  ）
（３）
（
２σ
√２πσ・成分は４６４個，成分は５７６個である．成分のを示したものである．この日の観測基線値は標本数回数，個数が多いのは，再観測が多かったためであが８個しかないため，ばらつきの大きい分布なのる．
か，それともばらつきの小さい分布に異常値が含ま図１を見ると，中央値は永年変化や季節変化を反れているのか，はっきりしない．
映して時間とともに変動しているが，一回ごとの観オミットのようなロバスト操作は，
が正規分布
測基線値の集団は中央値のまわりのおおよそ同じよに従うことを仮定し，分布に合わないデータを異常うな範囲に分布している．また，時折集団から外れと見なして除去することにより，残ったデータを正規分布と整合的にする操作と考えられる．しかし，
がどのような分布に従うかはっき
図２のように
りしない場合は，正規分布を仮定した異常値検出では精度が下がることが予測される．例えば，図２の分布の場合， ＝
３４．
２５，＝０．
１２なので，最も小さい観測基線値標本は  から の２．
５倍離れていることになるが，より正規分布に近い標本集団での の値０．
０７（第４章参照）を使うと４．
３倍となる．仮に，式２の を４とおくと，の値をどちらにするかによって，同じ標本でも異常かどうかの判断が分かれることになる．
図１から，観測ごとの中央値は時間変化するが，図２ ３月３１日に観測された H 成分の観測基線値の確率密各回の中央値と 個の観測基線値の差は時間によら度関数．

図３ ２０００年２月２５日〜１２月３１日に観測された H 成分の観測基線値のヒストグラム（左）と確率密度関数（右）．



地磁気絶対観測における基線値の異常値判定基準の定量化に向けてず同程度のばらつきをしていることが見て取れる．
そこで，回目の観測基線値を ，個数を 
，中央値としたとき，残差Δをすべての とを


＝ 
− 


について計算し，一つの標本集合としてみたときの確率分布（ヒストグラム），及び確率密度関数を調べた．図３に 成分の場合を示す．ヒストグラム
は０．
１の間隔で個数を計上している．
図３の標本数は４６４個あるので，図２に比べて分布が鮮明になった．０を中心とした釣鐘型の分布に少数の振幅の大きいサンプルが重ねあわさった形をしており，確率分布が式３で定義される確率密度関数でおおむね表現できていることから，おおよそ正規分布的であることがわかる．これは，観測エラー
の見積もりΔが，少数の異常な場合を除いて，観測者の違いや季節の変化によらず均質である，言い換えれば，同じ分布に従うランダム変数であることを示唆している．

図３から，残差Δ
の一年分の集合を使うことで標本数の不足を補い，残差分布のロバストなパラメータ推定を行う方法が考えられる．

５４）．残差が正規分布に従うとき，プロットの
Δ



は直線上に分布することが知られている
（藤井・シュルツ，１９９９）．異常値は，直線から大きく外れた値としてプロットされる．
図４において，３成分とも図の中央部は直線的で，順番の最も小さい部分と大きい部分で曲線的，あるいは断続的になっている．明らかに異常と思われるものだけを選びだすような保守的な検出を目指

４．解析
この章では，前章で説明した観測基線値データ中に含まれる異常値を検出するため，第２章のロバスト推定法を観測基線値データに適用する．前章で示されたように，一回ごとの観測基線値の集合では，標本数が少ないため分布のパラメーターの推定精度が下がり式２の異常値検出法がうまくいかない可能
性があるが，残差Δ
の一年分を使えば，標本数の不足を補って実効的な異常値検出ができると期待される．そこで，一年分の残差Δ
をを使って 求め，式２を変形して，

｜
Δ


｜


＞
 （４）を満たす残差中の異常値を検出することを試みる．
このとき，Δ
はすでに中央値が０になるよう平行移動された値であるから，左辺上部の中央値との差の演算は省略できる．
ま ず，の 値 を 計 算 す る．一 年 分 の Δ
のは，，，成分でそれぞれ，０．
０７，０．
０２２分，０．
０６であった．
次に，式２でしきい値を決めるパラメーター の値を決定する．このときに，異常値やそのしきい値のおおよその値を視覚的に捕らえるための簡便な手段として 





プロット（プロット）がある．プロットとは，で規格化した残差を小さい順に並べ，横軸に順番，縦軸にΔ
規格化した残差の値をプロットしたものである（図図４ H（上），D（中），Z（下）成分の QQ プロット．

６

伊藤信和・藤井郁子

図５

H（上），D（中），Z（下）成分の観測基線値と異常値．観測基線値を黒丸，検出された異常値を白丸で示す．



地磁気絶対観測における基線値の異常値判定基準の定量化に向けて

７

すとすると，プロットの連続的に並んだ標本かでは埋もれていた異常値の検出ができる，という２らずれているものが候補になる．・成分は残差 つの利点が得られる可能性が浮かび上がる．
の大きさがおよそ７，成分はおよそ５で連続的でなくなっている．そこで ・成分は ＝７，５．２ しきい値について
成分は ＝５をしきい値に設定した．これは，・
本報告では プロット図での連続性を考慮し・成分の残差振幅が，それぞれ０．
４９，０．
１１分，て，しきい値を ・成分は ＝７，成分は ＝０．
４２を越えるものを異常と判定するのと同値で５と設定したが，しきい値決定により客観性を持たある．
せるためにはさらに検討が必要である．
このようにして異常値の検出を行うと，成分は対象としているデータの性質にもよるが，藤井・８個，成分は１３個，成分は１２個の計３３個の残差シュルツ（１９９９）は，しきい値は経験的に３〜１２がが異常値と判定された．図５に，観測基線値と異常採 用 さ れ る と 述 べ て い る．一 方，（１９８９），徐値を示す．明らかに異常と思われる観測基線値デー（１９９８）は，そのしきい値は確率的には適切な意味タだけでなく，その日の分布状況では判断するのがを持たず任意に選ばれることから，ケース毎のロバ難しいデータも検出されていることが分かる．
スト推定の適用には限度があると主張している．しの値を小さくするほどしきい値が小さくなり， かし，第４章で述べたように，しきい値に対して異異常と判定される値が増える．ここで用いたしきい常値の影響を除いた残差の分布からおおよその意味値がどのような意味を持つか知るために，３成分そ付けを与えることも可能である．
れぞれの残差から検出された異常値を除いて 回目
本報告で扱った２０００年の観測基線値データに関しと全体の標準偏差σを求め 
を調の平均値μては，しきい値 ＝５〜７が一つの目安になると考−μ
べてみた．３３個のうち，成分７個，成分１３個，えられる．オミット基準として実用化するために
成分８個の計２８個の異常値が±３σの範囲を超えは，他の年のデータも調べて の値を確定する必要た．これにより，本報告で採用した異常値判定は， がある．
３σエディットルールと同程度かやや厳しい基準に基づくものであることが示唆される．
５．３ 標本数について
５．考察
５．１ 実際のオミット状況との比較
第４章での異常値検出結果と実際に行われたオミットとを比較して，本報告での異常値検出法と現行の方法との差を考察する．
本報告で検出された３３個の異常な観測基線値の実際の作業における取り扱い状況は，次の通りであった．
再観測をせずにオミットされたのは 成分１個（図４中の−１９．
５４）と 成分１個（図４中の−２１．
００），再観測の対象になりオミットされたのは１１個，オミットされなかったのは２０個，の計３３個である．
なお，本報告で異常と判定されなかった観測基線値の中で，成分の観測基線値データ（図４中の１．
２０９）が単独でオミットされているが，理由は不明である．
オミットされなかった場合を見てみると，その観測日には，しきい値に近い観測基線値データが複数存在している場合が多かった．そのため，一回の絶対観測の中では判定できず，残っていたものと考えられる．
以上により，本報告の異常値判定法を用いると，判定の自動化・客観化に加えて，（１）再観測をせずに異常値判定ができる，（２）１回ごとの絶対観測前節までに，ロバスト統計の適用により異常値を自動的に検出できることを示した．異常値を現行のオミットの対象とすれば，より真の値に近い基線値を求めることができると考えられるが，その際に問題となるのは標本数の減少であろう．第３章で見たように，１回ごとの観測基線値の数は正規分布を適用するには少ない状態にある．さらに標本数が減少すれば，最終的に得られる平均の基線値には大きなエラーバーがつくことになってしまう．オミット基準の策定に加えて，観測に必要な精度を確保するには最低何個の観測基線値が必要であるかを見積もり，オミットによりそれを下回った場合には再観測をするなどの規則の検討も必要であろう．
参考として，絶対観測で要求される精度を９５％の確率で達成するのに必要なサンプル数 を計算してみた．大きさ の標本で標準偏差が と分かっているとき，真の値 と推定値μとの最大誤差 μは
９５％の確率で１．
９６√より小さくなる．成分の基線値の許容できる最大誤差を０．
１０，はわからないので異常値を除いた標本標準偏差σ＝
０．
１６で代用すると，９５％信頼限界の条件を満たす は，
√

１．
９６
×０．
１６

０．
１０∴ ９．
８３（５）

８

伊藤信和・藤井郁子となる．従って，絶対観測の 成分が正規分布に従うとして，観測に要求される精度を９５％の確率で満たすには，少なめに見積もって１０個のサンプルが必要になる．精度を下げれば の値はもう少し小さくてもよい．
６．まとめ
鹿屋出張所の ，，成分の観測基線値１年分を使って，（１９８１）のロバスト統計の考え方で異常値を検出することを試みた．
異常値検出には，標本集団から中央値と を求め，各標本と中央値の差を で規格化することで，各標本が分布の中心からどのくらい離れているかを客観的に評価する．しかし，本報告で扱った
絶対観測では，１回分の観測の標本数が一般に８個程度と少なく，中央値と の計算精度が悪くなるので異常値検出の精度があがらない恐れがあった．そこで，観測基線値と各回の中央値との差を１年分集めたものを標本集団として，異常値計算を行った．
１年分の残差は ・成分は４６４個，成分は５７６個の標本で，０を中心としてほぼ正規分布し分布の端に少数の標本があった．，，成分のそれぞれの は０．
０７，０．
０２２分，０．
０６であった．
プロットにより，分布の連続性が途切れるところを異常値判定のしきい値と設定したところ，，成分で中央値からの差が の７倍，成分で５倍以上を異常とみなすこととなった．
このようにして検出された異常値は，，，成分でそれぞれ，８個，１３個，１２個あった．異常値を除いてから計算した残差の平均と標準偏差を使って，異常値と平均との差を標準偏差で規格化してみると，３と同程度かやや小さいしきい値を用いたのと同等であった．従来の異常値判定と比較すると，１３個はどちらでも異常と判定されていた．今回新たに異常とされた２０個について調べてみると，１回の絶対観測中に複数のしきい値に近い観測基線値が存在している場合が多かった．
このように１年分の観測基線値データを使って基準を設定して異常値の判定をしてみると，１日ごとでは埋もれていた異常値の検出も可能になり，ロバスト推定がオミット作業の効率化に役立つ可能性があることが示唆される．今後の課題としてもっと多くの事例を解析して，観測基線値に対するこの異常値検出法の客観性をより高めていくことが求められる．
参考文献
藤井郁子，シュルツ，
，地球磁場観測ネットワークデー
タの解析手法について（その１），研究会論文集，
９７
１
０４，１９９９．
     
３０８１９８１ ，非 線 形 フ ィ ル タ ー で デ ー タ を 洗 浄 す る，

，２００２．
徐 培亮，ロバスト推定の信頼区間と定数の重要性，地球惑星科学関連学会合同大会予稿集，９３，１９９８．
 
 

 １

３７
４

０

１９８９



地磁気絶対観測における基線値の異常値判定基準の定量化に向けて９
     
                ３１２００３３２００３７２００３                       ２５
２０００３１
２０００ 
５８５８

７２            ４６４  ５７６       ４６４５７６                            ７ ５    ８１２１３                     ２０００５７          