ノイズ除去方法とその方法を利用可能な生体情報計測装置および脳磁場計測装置
【課題】生体計測データから時間相関のあるノイズを除去することは難しい。
【解決手段】生体信号計測部１０は、生体から発生する脳磁場信号を計測する。ノイズ発生時刻取得部は、脳磁場信号を計測する際の心拍の発生時刻に関する情報を取得する。生体計測信号記憶部５０は、計測された脳磁場信号を心拍の発生時刻に関する情報とともに格納する。ノイズ信号抽出部２０は、脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出する。ノイズ除去処理部４０は、心拍の発生時刻に合わせて、抽出されたノイズ信号を脳磁場信号から減算することにより、ノイズ信号が除去された脳磁場信号を取得する。
【解決手段】生体信号計測部１０は、生体から発生する脳磁場信号を計測する。ノイズ発生時刻取得部は、脳磁場信号を計測する際の心拍の発生時刻に関する情報を取得する。生体計測信号記憶部５０は、計測された脳磁場信号を心拍の発生時刻に関する情報とともに格納する。ノイズ信号抽出部２０は、脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出する。ノイズ除去処理部４０は、心拍の発生時刻に合わせて、抽出されたノイズ信号を脳磁場信号から減算することにより、ノイズ信号が除去された脳磁場信号を取得する。
【発明の詳細な説明】
【技術分野】
【０００１】
この発明は、生体計測データからノイズを除去する方法、ならびにその方法を利用可能な生体情報計測装置および脳磁場計測装置に関する。
【背景技術】
【０００２】
生体計測データには必ず、取得したい信号とノイズ成分が混合されている。たとえば、人間の脳磁場を測定すると、計測データには、脳磁場の信号以外に、生体現象による心電、商用電源から発生するハムノイズ、およびガウスノイズが混在している。これらのノイズのうち、心電とハムノイズは時間相関のあるノイズであり、ガウスノイズは時間相関のないノイズである。生体計測データから所望の信号を得るためにはノイズ除去が不可欠である。
【０００３】
従来行われてきたノイズ除去の方法として、計測を繰り返して加算平均を取ることで、ノイズ成分を除去する方法や、フィルタリング処理を施す方法、独立成分分析（Independent component analysis、ＩＣＲ）によって計測対象の信号と独立なノイズ成分を分離する方法などがある。
【発明の開示】
【発明が解決しようとする課題】
【０００４】
急峻な成分を持つ心電や心磁などの生体現象ノイズは、加算平均やフィルタ処理を行っても十分に減衰しない。また、ハムノイズの除去にはバンドエリミネーションフィルタ（ＢＥＦ）が用いられるが、ＢＥＦは、商用電源の周波数の成分を全て減衰させるため、元の信号成分も減衰してしまう。また、ＢＥＦでは、商用電源の高調波成分を減衰させることができず、ノイズとして残ってしまう。
【０００５】
心電や心磁などの生体現象ノイズは、測定対象の脳磁場の信号とは独立であるから、独立成分分析をすれば原理的には分離することが可能であが、独立成分分析は、勾配法などにより繰り返し探索をする必要があり、分析には非常に長い時間がかかる。また、独立成分分析の適用条件が厳しいため、実際の測定が条件をすべて満たすことは難しく、実用的ではない。
【０００６】
本発明はこうした課題に鑑みてなされたものであり、その目的は、生体計測データから効率良くノイズを除去する方法、ならびにその方法を利用可能な生体情報計測装置および脳磁場計測装置を提供することにある。
【課題を解決するための手段】
【０００７】
上記課題を解決するために、本発明のある態様のノイズ除去方法は、生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとを含む。
【０００８】
「時間相関性のあるノイズ信号の発生時刻に関する情報」は、ノイズ信号の発生時刻のタイムスタンプ情報であってもよく、ノイズ信号の発生時刻（トリガ）を求めるための基礎となるデータであってもよい。たとえば、ノイズ信号の大まかな波形データがあれば、ノイズ信号の値が所定の閾値を超える立ち上がり時点を検出したり、ノイズ波形のピークポイントを検出することで、ノイズ信号の発生時刻（トリガ）を求めることができる。
【０００９】
本発明の別の態様は、生体情報計測装置である。この装置は、生体から発生する信号を計測する計測部と、前記計測部による生体計測信号の計測の際に発生する時間相関性のあるノイズ信号の発生時刻に関する情報を取得する取得部と、前記計測部により計測された生体計測信号を前記時間相関性のあるノイズ信号の発生時刻に関する情報とともに格納する記憶部と、前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するノイズ信号抽出部と、前記時間相関性のあるノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得するノイズ除去処理部とを含む。
【００１０】
本発明のさらに別の態様は、脳磁場計測装置である。この装置は、生体から発生する脳磁場信号を計測する計測部と、前記脳磁場信号の計測時の心拍の発生時刻に関する情報を取得する取得部と、前記計測部により計測された脳磁場信号を前記心拍の発生時刻に関する情報とともに格納する記憶部と、前記脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出するノイズ信号抽出部と、心拍の発生時刻に合わせて、抽出された前記ノイズ信号を前記脳磁場信号から減算することにより、前記ノイズ信号が除去された脳磁場信号を取得するノイズ除去処理部とを含む。
【００１１】
なお、以上の構成要素の任意の組合せ、本発明の表現を方法、装置、システム、コンピュータプログラム、データ構造、記録媒体などの間で変換したものもまた、本発明の態様として有効である。
【発明の効果】
【００１２】
本発明によれば、生体計測データからノイズを効率的に除去することができ、被験者の負担を大幅に軽減することができる。
【発明を実施するための最良の形態】
【００１３】
本発明の実施の形態の概要を説明する。実施の形態では、生体計測データから時間相関のあるノイズ除去を行う方法を提供する。生体計測データから、時間相関のあるノイズを抽出して除去する。時間相関のあるノイズは、ガウスノイズとは違って、加算平均だけでは除去されないことが多い。そこで、実施の形態では、時間相関のあるノイズをノイズ信号自身に同期させ、いったんノイズ信号を計測データから抽出し、抽出されたノイズ信号を計測データから減算する。
【００１４】
図６に時間相関のあるノイズの自己相関が示されている。左図は周期性のあるハムノイズの自己相関であり、１６ｍｓの周期があることがわかる。右図は心電波形の自己相関であり、約７００ｍｓ付近に１回目のピークがあり、ピーク値は１．０より小さくなっている。また、２回目のピークは約１４００ｍｓにあるが、ピーク値は１回目よりもさらに低減している。これは、心電が、ハムノイズのように正確な周期的ノイズではなく、時間的に揺らいでいるためである。
【００１５】
生体計測データには、計測対象の信号にノイズ成分が混合されている。図２に一例を示すように、生体信号には、時間相関のあるノイズとして、生体現象による心電と商用電源から発生するハムノイズが、時間相関のないノイズとしてガウスノイズが混合されている。
【００１６】
心電や心磁など急峻な成分をもつ生体現象ノイズはフィルタ処理を施しても十分に減衰しない。また、通常、ハムノイズ除去にはバンドエリミネーションフィルタ（ＢＥＦ）が用いられるが、ＢＥＦは商用電源の周波数を減衰させることができるが、その高調波を減衰させることはできない。図６８に示すように、ＢＥＦを用いると電源周波数である６０Ｈｚの成分すべてが減衰してしまい、また、※印で示す高調波が残ってしまう。
【００１７】
一般に周波数フィルタリングをかけると、ノイズ成分だけではなく、原信号の周波数成分も捨ててしまうことになる。図３５に示されるとおり、心磁ノイズは０．５Ｈｚ〜４０Ｈｚまでの周波数成分を含み、これは、計測対象の脳磁場（数Ｈｚ〜数十Ｈｚ）とスペクトルがオーバーラップしている。このため、周波数フィルタリングの手法では、心磁成分のみを選択的に除去することは原理的に不可能である。直観的には、心磁を除去すると、脳磁場も除去することになるからである。
【００１８】
そこで、本手法では、計測時に、時間相関のあるノイズの時刻情報（以下「トリガ」という）と計測データとを同時に記録する。計測データには、計測信号とトリガが含まれている。計測信号をノイズのトリガの周期の幅で区切り、区切られた計測信号をトリガのタイミングで同期させて、加算平均を行うと、時間相関のあるノイズのみが抽出される（図８の左）。計測信号を心拍の時刻情報に同期させて加算平均を取ると、心磁成分だけが強調されて抽出される。ここで、心拍は時間的にゆらぎがあることに配慮して、加算平均時間の前後に余裕をもたせ、加算平均を取る時間の範囲を一部重複させている。加算平均を取る時間の範囲を重複させる幅は、心拍トリガが時間軸方向にゆらぐ量に応じて決める。このように加算平均を取る時間の範囲に余裕をもたせておくことで、心拍トリガのタイミングが多少前後にずれても、心拍トリガの周期で区切られた計測信号の間で心拍トリガが一致するように時間軸方向にずらして加算平均を取ることができる。
【００１９】
心電信号が急峻に変化するタイミングをトリガとして検出するが、トリガの検出方法としては、心電信号の値が閾値を超える立ち上がりポイントを検出すればよく、実施の形態ではこの方法が用いられる。トリガの別の検出方法は、心電信号のピークを検出することである。いずれにしても大まかな心電波形さえ得られれば、その心電波形から急峻な立ち上がりポイントをトリガとして検出することで心拍の時刻情報を取得することができる。
【００２０】
次に、加算平均によって抽出したノイズ成分を計測信号から除去する（図８の右）。具体的には、同時記録した心拍の時刻（心拍トリガ）に合わせて抽出されたノイズを元データから減算する。これにより計測された信号から心磁ノイズが除去される。このとき、加算平均前後に時間的余裕を持たせているので、心磁ノイズのような時間的にゆらぎを有するノイズにも対応できる。
【００２１】
また、本手法をハムノイズ除去に適用するとＢＥＦとは異なり、図６９に示すように、電源周波数の信号成分を減衰させることなく、高調波成分を含めてハムノイズだけを除去することができる。
【００２２】
このように、本手法は、ハムノイズ等の周期的なノイズの除去、生体現象による時間的にゆらぎのあるノイズ（たとえば、心拍ノイズ）の除去の両方に対応することができる。
【００２３】
本手法によって時間相関のあるノイズを生体計測データからあらかじめ除去することができ、残りのガウスノイズなどは加算平均によって除去すればよいから、通常のノイズ除去に必要な加算平均回数を半分から４分の１程度に減らすことができる。そのため、被験者の負担を大きく軽減することができる。たとえば、被験者に音を１秒に１回聴かせて脳磁を計測する実験を、従来は１０００回くらい繰り返して加算平均を取っていたため、実験には２０分近くかかっていた。これが少なくとも半分から４分の１程度の加算平均回数で済むため、実験は５〜１０分程度の時間に短縮される。
【００２４】
本手法は、脳波計測やＭＥＧ計測におけるハムノイズ除去や心磁ノイズ除去に応用できることを実施の形態で詳細に説明する。しかし、これは一実施例に過ぎず、本発明のノイズ除去方法は、脳波やＭＥＧ計測に限らず、一般に生体計測データから時間相関性のあるノイズを除去することを目的とするあらゆる技術分野に適用することができる。
【００２５】
図１は、本発明の実施の形態に係る脳磁場計測装置１００の構成図である。脳磁場計測装置１００は、後述のＳＱＵＩＤセンサを用いた脳磁図（ＭＥＧ）測定装置であるが、ここでは、脳磁場を測定するための一般的な構成については省略し、ノイズ除去に関する構成について説明する。
【００２６】
生体信号計測部１０は、生体から発生する脳磁場信号を計測して、生体計測信号記憶部５０に記録する。ノイズ発生時刻取得部３０は、同一の生体から測定される心電波形をもとに心拍の発生時刻に関する情報を取得して生体計測信号記憶部５０に記録する。生体信号計測部１０による脳磁場信号の記録と、ノイズ発生時刻取得部３０による心拍の発生時刻の記録は同時に行われる。
【００２７】
ノイズ信号抽出部２０は、生体計測信号記憶部５０に記録された脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出し、抽出されたノイズ信号をノイズ信号記憶部６０に記憶させる。
【００２８】
ノイズ除去処理部４０は、心拍の発生時刻に合わせて、ノイズ信号記憶部６０に記憶されたノイズ信号を生体計測信号記憶部５０に記録された脳磁場信号から減算することにより、ノイズ信号が除去された脳磁場信号を取得し、ノイズ信号除去後の脳磁場信号を生体計測信号記憶部５０に書き戻す。
【００２９】
計測された脳磁場信号には、心拍が原因となって発生する心磁成分がノイズとして混在する。本手法では、心磁信号を直接測定することなく、心電波形が得られる心拍の発生時刻情報だけを用いて、計測された脳磁場信号を心拍の時刻に同期させて加算平均を取ることにより、計測された脳磁場信号から心磁信号を抽出することができる。抽出された心磁信号を計測された脳磁場信号から差し引けば、心磁信号を除去して脳磁場信号だけを取り出すことができる。
【００３０】
計測対象の生体信号は、脳磁場信号に限られない。計測対象は、脳波信号であってもよく、その他の生体から発生する任意の信号が計測対象となりうる。また、計測対象の生体信号にノイズとして混在する信号は、心電信号に限らず、時間相関のある任意の信号がノイズとなりうる。たとえば、計測装置の交流電源から発生する周期性のあるハムノイズがノイズ信号として生体計測信号に混在する場合にも同様のノイズ除去手法を適用できる。すなわち、ハムノイズ信号の発生時刻情報を取得し、発生時刻に同期させて生体計測信号を加算平均することでハムノイズ信号を抽出することができる。抽出されたハムノイズ信号を生体計測信号から差し引けば、ハムノイズのない生体計測信号が得られる。なお、交流電源から発生するハムノイズの発生時刻は、交流がマイナス方向からプラス方向にゼロを切るタイミングを検出することで得られる。この発生時刻情報の記録を生体信号の計測と同時に行えばよい。
【００３１】
以下、本発明の実施の形態を発明者の「心磁成分が脳磁場計測データに及ぼす影響」に関する研究内容に即して詳細に説明する。脳磁場計測装置１００の各構成の動作の詳細は、下記の研究内容を参照して理解される。
【００３２】
本研究は、ＭＥＧ計測において生体現象ノイズの心磁成分に注目した。心磁成分は、カットオフ周波数３［Ｈｚ］程度のフィルタ処理だけでは減衰が困難である事と、注目しているＥＣＤの位置によって、その影響が異なることが分かった。そこで、加算平均を行い心磁成分の抽出を行い、これを除去するシステムの開発を行った。聴覚と視覚のＭＥＧ実験を行い、本システムの評価を行った結果、良好な性能を有することが示された。
【００３３】
以下の説明の章立てを示す。
１章では、ＭＥＧ実験中に計測された様々な時間相関のあるノイズの概観を解説し、心磁成分除去について問題を設定した。また、本手法で提案する心磁成分の除去手法について解説する。
【００３４】
２章では、ＭＥＧで使用できる心拍計測システムについて解説する。本研究で行う心磁成分の除去手法は、心拍をＭＥＧ計測中に同時に記録しなくてはならない。そこで、ＭＥＧ計測中に問題の起こらない心拍計測装置の開発と既存の計測装置との評価を行った。
【００３５】
３章では、ＭＥＧによる心磁成分の抽出を行った結果について解説する。抽出した心磁成分の特徴を明らかにし、従来行われているフィルタ処理等では減衰が困難だということを示した。また、人工的なデータを用いて、推定部位によって心磁の影響が異なることを示唆した。
【００３６】
４章では、既存の脳磁図解析ソフトウェアでは出来ない本手法による心磁除去を実現するプログラムについて解説する。
【００３７】
５章では、１次聴覚野、１次視覚野の脳磁場応答が計測できる実験を設計し本手法の評価を行った。その結果、心磁成分が加算平均により、基線の揺らぎの一要因となっていることを示し、本手法の有効性を示すことが出来た。
【００３８】
１．序論
ＭＥＧ計測はノイズとの戦いである[1]。ＭＥＧ（Magneto Encephalo Graphy：脳磁図）は、大脳皮質に生じる微弱な電流から生じる磁場を多数のＳＱＵＩＤ（Superconducting Quantum Interference Devices：超伝導量子干渉素子）センサを用いて計測する装置である。誘発脳磁場は大変微弱である為に、磁気シールドなどで外部の磁気ノイズを減衰させて計測を行う。しかし、この磁気減衰効果は、生体から発生する磁気ノイズに対しては働かない。その為に、計測される脳磁場データには、目的以外の信号が多数含まれる（図２）。通常は計測された脳磁場データに対して、加算平均や信号処理などの手法を用いて、Ｓ／Ｎ比（Signal/Noise Ratio）を向上させる。
【００３９】
計測されるノイズには、大きく分けて生体現象が原因のものと計測環境が原因のものの二つがある（図３、図４）。
【００４０】
図２は、ＭＥＧで計測される信号を説明する図である。計測される信号は、ノイズ成分、心磁成分（生体現象によるノイズ）、ハムノイズ（環境に起因するノイズ）等が、線形に混合された磁場信号が計測されている。
【００４１】
図３は、ＭＥＧで計測された生体現象によるノイズを説明する図である。１６０個のＳＱＵＩＤセンサから計測された磁場信号を重ね合わせた磁場波形と等磁場線図。等磁場線図の緑は磁場の吸い込み（Sink）、赤は磁場の吹出（Source）を示す。等磁場線図は黄線の時刻のものである。ａ）の眼球運動はＳＱＵＩＤセンサに近い為に大きな磁場変化がみられる。ｅ）の呼吸は、非常にゆっくりとした変化である。また、等磁場線図も他とは大きく違う。
【００４２】
図４は、ＭＥＧで計測された環境が起因となっているノイズを説明する図である。ＳＱＵＩＤセンサから計測された磁場信号とパワースペクトル分布を示す。６０［Ｈｚ］、１２０［Ｈｚ］、１８０［Ｈｚ］というハムノイズの６０［Ｈｚ］の倍数の周波数においてピークが見られる。このピークは、６０［Ｈｚ］の交流電源から発生していると考えられる。
【００４３】
ノイズには、時間に対して相関を持たないものと時間に対して相関を持つものがある（図５、図６）。誘発磁場以外の信号が、時間に対してランダムに変動すれば（時間と無相関な信号であれば）、Ｎ回の加算平均で１／Ｎの振幅に減衰する。これに対して、時間に対して相関のあるノイズ（時間相関のあるノイズ）信号はランダムに変動しないので、加算平均と同期してしまい、いつまでも相殺されずに残る可能性がある。また、生体現象ノイズ信号は通常脳磁場の１０倍程度の大きな信号である事が多いので、加算平均回数が１０回程度でも脳磁場と同程度の影響を結果に及ぼす。そのため、時間相関のあるノイズは、できるだけ加算平均から除くようにしなければいけない。一方時間相関があることを利用して、ノイズ成分を抽出できる可能性がある。
【００４４】
図５は、ランダムなノイズの自己相関を説明する図である。計測によって生じるノイズは、時間に対して相関を持たない。その為に、加算平均を行うと減衰する。
【００４５】
図６は、時間相関のあるノイズを説明する図である。左図は周期性のあるハムノイズの自己相関である。赤丸で示すようにピークが周期的に一定の間隔で生じている。右図は、心電の自己相関を表す。Ｒ−Ｒ間隔の部分を赤線で示す。心臓の拍動は一定の間隔で生じるが、時間的に揺らぎを生じる。
【００４６】
ＭＥＧ実験を行うとき、図７のように被験者の動きを抑制する。頭部が動くことによって生じる生体現象ノイズを抑制する。また、注視してもらうことで眼球運動が起こりづらくなるような実験タスクの工夫を行う。それでも、呼吸や心臓から、生体現象ノイズは生じてしまう。しかし、呼吸や眼球運動の生体現象のノイズは、通常計測される磁場データより数倍波形が大きく変化するので、加算平均時に含めないようにすることが可能である。また、電源から生じるハムノイズは計測時のバンドエリミネーションフィルタ（ＢＥＦ）、計測後のフィルタ処理により、大幅に減衰させることが可能である。しかし、心磁成分は、一定の間隔で生じ続け、除外することは困難である。
【００４７】
図７は、MEG実験中の風景を説明する図である。デュワーと被験者の頭部との間を、弾性がある素材で埋めている。これにより、頭部の横の動きは抑制される。また、野村が製作した顎乗せ台[4]を用いて、頭部の上下の動きを抑制している。顎乗せ台を使用する前のマーカーコイルずれは最大で９．７２［ｍｍ］であり、顎乗せ台を使用した後のマーカーコイルのずれは最大で０．４０［ｍｍ］になったという報告がある。
【００４８】
ＰＣＡやＩＣＡといった信号処理の手法による、ノイズ除去も提案されている[5][6]。独立成分分析（ＩＣＡ）は複数の観測信号から統計的な独立性を元に、源信号に分離する手法である。ノイズと誘発磁場は無相関である為に、互いに分離が可能である。しかし、ＭＥＧは多数のＳＱＵＩＤセンサを用いて計測するために、多次元のデータになり、解析には多くの時間がかかる[7]。
【００４９】
本研究の目的は、生体現象ノイズの心磁成分に注目し、心磁成分を除去するシステムの開発と等価電流双極子推定への結果に及ぼす影響を調べることで、このシステムの性能評価を行うことである。
【００５０】
本研究で行う心磁成分の除去手法は、心拍をトリガにし、加算平均を用いて心磁成分の抽出を行い、除去を行うという手法である。心磁成分は、一定のリズムで発生する為に、心拍を同時に記録すれば、心磁成分のみを除去するのに有効だと考えた。
【００５１】
１．１．加算平均による減算手法
本研究で提案する心磁成分の除去手法は、心臓の動きを、脳磁場計測と同時に記録する。これをトリガに用いて加算平均を行い、心臓の拍動によって生じる心磁成分を抽出し、オフラインで計測データから減算を行うという手法である（図８）。
【００５２】
一般に運動負荷を与えると心電波形のＱＲＳ時間、ＱＴ間隔が短縮する[8]。複写や計算課題などの負荷を与えると、Ｒ−Ｒ間隔が変動し精神的負荷を反映するといわれている[9]。本研究では、実験中ほぼ安静である為に抽出される心磁波形がほぼ相似になると考えた。Ｒ−Ｒ間隔の変動に対しては、１周期分の心磁波形を求め、心拍トリガから除去する事により対応が出来る。
【００５３】
図８は、本研究が提案する心磁成分を除去する手法を説明する図である。まず心拍の時刻情報を元に、ＳＱＵＩＤで計測された信号に対し、加算平均を行う。加算平均された信号は、心磁に由来する磁場成分であると考えられる。その後、加算平均に用いた心拍の時刻を用いて、抽出した心磁成分で除去を行う。
【００５４】
本手法を実現する為に、ＭＥＧ計測に影響を及ぼさない心拍計測装置と、本手法を行うソフトウェアが必要になる。
【００５５】
１．２．基線補正
加算平均を行ったデータに対し基線補正をかける。基線補正は、なんらかの持続的な神経活動や基線の揺らぎの除去を行う。基線補正が行う操作を図９に示す。
【００５６】
図９は、基線補正を説明する図である。指定した時間区間の各ＣＨの平均値を求めて、１サンプルごと元データを減算する。刺激トリガの前から１００［ｍｓ］、全区間で行われる。この操作は、直接計測データを変化させる。
【００５７】
直接計測データを変化させる為に、平均を求める区間によって、本来得たい信号を消してしまう可能性がある。また、どの部分を使えばよいという基準が明確に決められていない。
【００５８】
１．３．心電と心磁について
心臓は心筋を収縮させることによって、血液を全身に循環させるポンプの機能をもつ臓器である。心臓には、左右に隔壁と弁膜で仕切られた心室と心房が備わっている。それ以外にも、心房と心室を構成する心筋を収縮コントロールするための電気的指令を出したり、伝えたりする組織を供えている（図１０）。
【００５９】
図１０は、心臓の各部位の名前と特徴を説明する図である。洞結節は心筋の収縮の一定のリズムを作り出す部位である。左右房の心筋の脱分極による伝播は三尖弁と僧房弁でとまり、左右室への伝播は房室結節のみに伝わる。心室の刺激伝導系はヒス束から左右室の心筋まではプルキンエ繊維抹消である。この間は非常に早く伝わる。
【００６０】
心筋細胞は、静止時に陰性に保たれ、興奮によって脱分極する。この際、Ｎａ＋が細胞内に流入する。細胞はこの流入したＮａ＋を積極的に排出する。３個のＮａ＋が排出されると、２個のＫ＋が流入する。この時、外向きに電流が生じる。心筋が収縮する時Ｃａ２＋が必要になる。心筋の脱分極直後は、Ｎａ＋が多いのでＮａ＋を排出し、Ｃａ２＋を取り込む（外向き電流が生じる）。Ｃａ２＋が多くなると、逆にＣａ２＋を排出し、Ｎａ＋を取り込む（内向き電流が生じる）。心筋が収縮後、次の脱分極に備えて、再分極が起こる。この時、細胞内のＫ＋を排出し、元の静止膜電位に戻ろうとする[12]。
【００６１】
このように、心筋細胞はＮａ＋により脱分極、Ｃａ２＋により収縮、Ｋ＋により再分極がおこわなれる。心筋の興奮伝導（脱分極の伝播）を説明する。静止時、細胞内は陰性に保たれているが、隣り合う心筋細胞により脱分極が伝播する。この時、心筋の細胞内は陽性になる。興奮伝導の波面に電流双極子を生じる[11]。
【００６２】
体表面に電位差が生じる仕組みを説明する。心筋の興奮、収縮に伴う電流があると抵抗性を持つ周囲の臓器に帰還電流が生じる。その為に、任意の点で電位差が計測される（参考文献[10]の図２）。
【００６３】
心電波形と心筋細胞の収縮タイミングを説明する。計測される心電波形は、心筋の収縮タイミングを表している。参考文献[10]の図３に示されるように、心電図は心臓の電気的な状態を示し、Ｐ波、Ｑ波、Ｒ波、Ｓ波、Ｔ波、Ｕ波の順に現れる。Ｐ−Ｑの間に、左右房の収縮が起こる。ＱＲＳの間には、ヒス束からプルキンエ繊維抹消への興奮伝導が起こる。つまり、心室筋の脱分極のはじまりと終わりを表す。Ｓ−Ｔは、次の収縮の為に再分極している。
【００６４】
心磁の電流源は、興奮伝導に伴う電流双極子だと考えられる。心筋の興奮の伝播が心磁となって計測されると考えられる。
【００６５】
２．ＭＥＧで使用できる心拍計測システムについて
心臓は、心筋の収縮によって全身へ血液を運ぶ臓器である。この変化は、血管内の容積変化（脈波）と心筋の筋収縮によって生じる体表面の電位変化（心電図）によって計測できる。本研究の心磁除去の為に、心拍と同期した信号を、ＭＥＧ計測中に同時に記録しなくてはいけない。そこで、ＭＥＧ計測と同時に計測しても、ノイズや実験タスクに問題の発生しない心拍計測装置について検討を行った。
【００６６】
２．１．心拍を計測するシステムの評価
ＭＥＧ実験では、体性感覚用の刺激装置に電気刺激装置がある。また、実験タスクによっては、被験者にボタン押し動作を行う事がある。このような刺激や動作の時に、本件急で開発する心拍を計測するシステムが影響をうけないようにしなければいけない（図１１）。
【００６７】
図１１は、実験タスクの例を説明する図である。体性感覚用の刺激装置に電気刺激装を用いるものがある。実験タスクによっては、ボタン押し動作を行うことがある。これらの実験中に、心拍計測するシステムに影響を与えないようにしなければならない。
【００６８】
２．２．脈波を用いた心拍の計測
脈波とは、心拍動によって生ずる動脈系圧波動の伝達である[13]。この伝達を計測する方法に、空気動圧センサを用いて変動を計測する方法とヘモグロビンの赤外光における吸収量が違うことを利用する方法がある。
【００６９】
空気動圧センサを用いた脈波を計測する装置に（株）エム・アイ・ラボ製の脈波観測装置（図１２）を用いた。赤外線を用いた脈波観測装置は、赤外線ＬＥＤと光センサ（フォトダイオード）で構成した心拍計測装置を製作した（図１３は回路図、図１４は製作した装置）。
【００７０】
図１２は、脈波観測装置（空気圧の変化で捉える）を説明する図である。脈波観測装置は空気袋のわずかな空気圧の変動を捉え、電圧に変換し出力している。使用時は、空気袋を被験者の指にサージカルテープで密着させて使用する。このときチューブが空気を伝わる遅延と心臓から指までの血流の遅延を考慮しなくてはいけない。
【００７１】
図１３は、脈波観測装置（赤外線を使用した装置）の回路図である。赤外線LEDとフォトダイオードの間に指を挟むことにより、血流の変化を計測することができる。
【００７２】
図１４は、製作した脈波観測装置を説明する図である。赤外線は、人の目で見ることが出来ない。左図で、赤外ＬＥＤ部分が光っている。これは、デジタルカメラの光センサが赤外光に感受性をもつため記録される。右は、脈波の測定風景を示す。赤外線LEDとフォトダイオードを隠すように指を置いて計測を行う。
【００７３】
２．３．心電計の製作
心筋の細胞が興奮することによって生じる電場は、ほぼ、瞬間的に周辺の臓器を経由して体表面まで到達する。心電計は、この電場変化を捉える装置である。臨床目的の心電図は、標準１２誘導（肢誘導と胸部誘導）と呼ばれる電極の貼り付けを行う[15]（図１５）。本研究では、心電を臨床目的に使用するわけではないので、単純に心筋の収縮によって生じるＱＲＳの電場変化を捉えて増幅すればよい（図１６）。製作した心電計の回路図と装置概観を図１７、図１８に示す。
【００７４】
図１５は、標準１２誘導の電極貼り付け位置を説明する図である。肢誘導と胸部誘導の電極の取り付け位置を示す。この電極貼り付け位置は、臨床目的で使用されている。
【００７５】
図１６は、典型的な心電波形を説明する図である。心筋の細胞が、脱分極すると、心筋は収縮をする。Ｐ波は、心房の収縮を表す。ＱＲＳ波は、心室の脱分極の始まりから終わり、つまり心室の収縮を表し、Ｔ波で心筋の収縮を終えることを表している。詳細は１．３で解説している。
【００７６】
図１７は、製作した心電計の回路図である。ＡＤ６２０は計装アンプである。高い入力インピーダンスを持った差動増幅専用のオペアンプである。ＶＲで、ゲインを１〜１００００倍まで増幅できる。ＲＣで構成された０．５［Ｈｚ］のハイパスフィルタと４０［Ｈｚ］のローパスフィルタをかけた後、非反転増幅し出力している。
【００７７】
図１８は、製作した心電計を説明する図である。製作した心電計は、シールドルーム外に設置する。ケーブルをシールドルーム内まで伸ばし、被験者に電極を取り付ける。
【００７８】
被験者の体につける電極は、日本光電製ディスポ電極ビトロードF-150Sを用いた。電極は、実験前に磁石を用いて、磁性体ではないことを確認している。市販の製品は、ディスポ電極と心電計をつなぐ電極リードにステンレスが用いられており、ＭＥＧ実験に使用できないことが分かった。その為、磁性体を用いていない電極リードを自作した。導線は、銅を用い、ビトロードとの接続部分には、電子工作用に用いられるゴム足を加工して用いた（図１９）。
【００７９】
図１９は、製作した電極リード（左）と日本光電製の電極リード（右）を説明する図である。自作した電極リードは銅線部分がディスポ電極に接触するように出来ている。日本光電製の電極リード（使い捨て電極リード線BC−IXX：XXは色によって型番が異なる）はディスポ電極をバネの構造で密着させている。
【００８０】
計測した心電データに対して、デジタルフィルタをかける。このフィルタ処理は、計測中の被験者の動きによる基線の揺らぎ（〜２［Ｈｚ］）や、ハムノイズ（６０［Ｈｚ］の成分）を除去する為のものである。
【００８１】
２．４．脈波観測装置と製作した心電計の心拍計測結果
図１２と図１４の脈波観測装置を用いて、心拍を計測した結果を示す（図２０、図２１）。この計測の被験者は、共に、２秒間に３回ほど、振幅のピークが見られた。安静時の心拍の間隔が、正常な間隔が１秒間に０．８〜１．７回[14]なので、正常に心拍が計測されているものと考えられる。
【００８２】
図２０は、空気動圧センサを用いた脈波計測の結果を説明する図である。指を動かすと、大きく基線がゆらいだ。また、空気袋を触ると大きな変化になってしまい、脈波が正常に観測されない。
【００８３】
図２１は、光センサを用いた脈波計測の結果を説明する図である。基本的に指を光センサと赤外線LEDとの間に置き続けなければ計測が出来ない。その為に、ＭＥＧ実験中に、ボタンを押すという動作は出来ない。
【００８４】
計測した心電波形に、３−５０［Ｈｚ］のＦＩＲ型のバンドパスフィルタをかけた。これにより、被験者が動くことによって生じる心電波形の基線の揺らぎとハムノイズを除去できる（図２２）。この時、位相ひずみが起きないようにフィルタは直線位相特性を持ったものを使用した[16]。
【００８５】
図２２は、フィルタをかける前の心電波形とバンドパスフィルタ後の心電波形を説明する図である。細かく一定の周期で現れているのがハムノイズである。交流電源の周波数が６０［Ｈｚ］なので、その帯域を減衰するように３−５０［Ｈｚ］のバンドパスフィルタをかけた。
【００８６】
製作した心電計は２点間の電位差を求めて増幅することが出来る（予備実験の計測は、ノイズ低減の為に不感電極を心電計のＧＮＤに接続して計測を行っている。）。よって、どの部位に電極を貼り付けるかによって、波形は異なる。そこで、電極を頂点とみなし、三角形に心臓が収まるように電極をつけて計測を行った（図２３、図２４）。
【００８７】
図２３は、胸部誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。赤がプラス入力、青がマイナス入力、黒が不感電極（心電計ＧＮＤ）に接続している。
【００８８】
図２４は、肢誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。赤がプラス入力、青がマイナス入力、黒が不感電極（心電計ＧＮＤ）に接続している。全体的に、胸部誘導より、振幅が低くなった。これは、心臓からの距離が離れた為に、臓器などが抵抗成分となった為だと考えられる。
【００８９】
ＭＥＧ実験で使用する電気刺激装置とボタン押し動作をした時の結果を図２５に示す。心電波形に変化は見られなかった。
【００９０】
図２５は、心電計測中に電気刺激、ボタン押し動作をしたときの心電波形を説明する図である。赤い線は、それぞれの刺激の始まりのタイミングを表している。電気刺激は０．５［ｍＡ］を手首に出力している。ボタン押し動作は、人差し指でボタンを押してもらった。共に心電波形に変化は見られなかった。
【００９１】
製作した脈波観測装置（図１４）と心電計を用いて同時に心拍を計測した結果を図２６に示す。心電のＲ波と脈波のピーク値はおよそ２００［ｍｓ］の遅延が生じていた。この遅延は、被験者によって心臓から脈波の測定部位までの距離が異なる為に、変化する。
【００９２】
図２６は、光センサを用いた脈波観測装置と製作した心電計の心拍計測を説明する図である。脈波は心電に比べて遅延が生じていた。これは、血流の変化を捉えている為である。
【００９３】
２．５．考察
本研究では、MEG計測と心拍を同時に記録し、加算平均することにより心磁成分のみを抽出する。その為に、実験タスク中でも安定して心拍を計測する装置が必要になる。
脈波観測装置は、空気動圧センサを用いるもの、赤外線を用いるもの、それぞれに問題がある。実験タスクによって、ボタンを押すという動作をすることがある。図２０から、空気動圧センサを用いた脈波観測装置は、使用できない。光センサを使用する脈波観測装置は、センサをシールドルーム内に持ち込むことは出来ない。光センサの電流の変化が、ＳＱＵＩＤセンサに混入してしまうからである。また、光ファイバを用いることも考えたが、赤外線の波長の減衰が大きく、脈波成分が計測出来なかった。
【００９４】
製作した心電計は、電気刺激、ボタン押し動作を行っても、安定して心拍が計測できた。また、脈波とは違い、電場の変化を計測しているので、遅延が生じない。そこで、本研究では、心拍を計測する装置に心電計を用いて、心磁成分の抽出を行うことにした。
【００９５】
３．ＭＥＧによる心磁成分の抽出
生体から発生する磁気を用いて機能を計測する装置は多数存在する。脳磁図は文字通り、脳機能をターゲットに作られた計測装置である。それ以外にも、心臓の機能を調べる心磁図用の計測装置がある。本研究においては、ＭＥＧ計測装置を用いた心磁図の報告は無いので、心磁がどのようにＳＱＵＩＤセンサに計測されるか調べた。
【００９６】
参考文献[17]に示される心磁図計測装置を用いる。心臓の真上に、ＳＱＵＩＤセンサアレイを配置し、計測を行う。ＭＥＧ計測装置は、頭部を覆うようにＳＱＵＩＤセンサを配置している。
【００９７】
心磁図測定装置を用いて計測された心磁成分は、参考文献[17]に示されるように、ＳＱＵＩＤセンサ８ｃｈ分を重ね合わせた結果である。
【００９８】
３．１．実験方法
被験者に電極を図２７のように貼り付けた。また、心電計の出力をＭＥＧ計測装置のトリガチャンネルに接続し、心電を同時に記録した。被験者は２２−２３歳男性４名で行った。
【００９９】
図２７は、電極の貼り付け位置を説明する図である。心臓が電極を頂点に、三角形の中に捉えられるように、反転入力、非反転入力、不感電極を貼り付けた。
【０１００】
心磁成分の抽出には、同時に記録した心電図のＲ波部分を、トリガに加算平均を行った。解析に用いたソフトは、脳磁図用の解析ソフトウェアMEG Laboratory（横河電機、イーグルテクノロジー製）である。加算平均は２００回以上行った。通常の脳磁場の解析に用いられるフィルタ処理、基線補正処理は行わなかった。
【０１０１】
３．２．ＭＥＧ計測装置について
脳磁場計測と心磁成分の計測には全頭型１６０チャンネル脳磁計（イーグルテクノロジー製）を使用した。ＳＱＵＩＤセンサはGradio-meterでBaselineは５０．０［ｍｍ］であった。ノイズとシグナルのピーク値を用いて計算したＳ／Ｎ比は０．０３［ｄＢ］であった。被験者は磁気遮断室（ＭＳＲ：Magnetically Shielded Room）内に入ってもらう。脳磁場計測には、顎乗せ台を用い、デュワーと頭部の間を図７のように埋めて計測を行った。
【０１０２】
図２８は、全頭型１６０チャンネル脳磁計を説明する図である。被験者は椅子に座り、デュワー内に頭部を入れてもらい脳磁場計測を行う。
【０１０３】
全ての計測時に共通なパラメータを表１に示す。
（表１）
ＭＥＧ計測のパラメータ
サンプリング周波数［Ｈｚ］ １０００
Low Pass Filterのカットオフ周波数［Ｈｚ］ ２００
Band Elimination Filterの中心周波数［Ｈｚ］ ６０
【０１０４】
３．３．ＭＥＧによる心磁成分の抽出結果
心電波形をトリガに、加算平均を行った結果を示す。図２９は加算平均によって得られた波形と、図３０はSQUIDチャネルの波形変化を示す。
【０１０５】
図２９は、加算平均によって抽出した心磁波形と等磁場線図を説明する図である。等磁場線図の緑は磁場の吸い込み（Sink）、赤は磁場の吹出（Source）を示す。変化の激しい時刻（赤線、青線）の等磁場線図を示す。
【０１０６】
図３０は、各ＳＱＵＩＤセンサに含まれる心磁成分を説明する図である。各ＳＱＵＩＤセンサを大きく左右に分けると、左側で最も大きく変化するチャネルは１０８ｃｈ、右側で最も大きく変化するチャネルは７１ｃｈであった。
【０１０７】
心電をトリガに、図２９のような波形が得られたことから、心拍による誘発磁場だと考えられる。同時に記録した心電波形と図３０から分かった変化が著しいチャネルの計測磁場データとの比較を行ったものを図３１に示す。
【０１０８】
図３１は、心電波形と計測磁場データとの比較を示す図である。心電波形のＱＲＳ波にかけて、７１ｃｈ、１０８ｃｈに同期した変化が見られる。７１ｃｈから１０８ｃｈを減算するとより分かりやすい。
【０１０９】
図３２は、心電波形と計測した心磁図の比較を示す図である。Ｐ波、ＱＲＳ波、Ｔ波に対応した心磁成分が計測されている。ＱＲＳ波に相当する部分が大きく計測された。
【０１１０】
被験者３名の心磁成分抽出結果を図３３に示す。全被験者でＰ、ＱＲＳ、Ｔに相当する心磁波形が確認できた。
【０１１１】
図３３は、被験者３人の心磁成分抽出結果を説明する図である。等磁場線図は図２９の青線と同じ心磁成分の位置を示す。被験者によって、ＱＲＳ波に相当する心磁成分の波形が異なっていた。また、振幅のピーク値も被験者によって異なった。被験者４人計測して、振幅の絶対値を取った所、最大で８２０［ｆＴ］、最小で３３３［ｆＴ］、平均して４８５±２２９［ｆＴ］の振幅が計測された。
【０１１２】
心磁波形が時間的に変化するか調べた（図３４）。心磁成分の周波数成分を図３５に示す。
【０１１３】
図３４は、心磁波形の相似を説明する図である。図のＡ、Ｂ、Ｃは、それぞれ時間の区間を表す。約２００［ｓ］で３００回の加算平均を行った。心磁波形にほとんど変化は無かった。
【０１１４】
図３５は、心磁波形のパワースペクトル分布（１０８ｃｈのみ）を説明する図である。ＱＲＳの心磁成分が急峻な変化である為に、広い周波数範囲に成分が含まれている。図によるとほぼ１〜５０［Ｈｚ］までの成分が含まれている。１〜５０［Ｈｚ］の成分は、脳磁場解析にも含まれる帯域である。
【０１１５】
３．４．心磁成分が脳磁場計測データに与える影響
抽出した心磁成分に対し、等価電流双極子推定をおこなった結果を図３６に示す。
【０１１６】
図３６は、抽出成分に対して等価電流双極子推定をおこなった結果を説明する図である。推定位置の色とＧＯＦ、Intensityがそれぞれ対応している。心磁波形の振幅が最も高い時刻に、４２０［ｎＡ・ｍ］ものIntensityが計測された。推定位置は、脳中央に推定された。推定アルゴリズムはSarvas法を用いた。
【０１１７】
心磁成分が、誘発磁場応答に、どのように影響を与えているか調べる為に図３７のような磁場データを人工的に作成した。作成にはMEG Laboratoryを用いた。手動で右の１次聴覚野近傍と１次視覚野近傍に等価電流双極子を置いた（図３８）。それぞれの等価電流双極子が作る磁場データに対し、心磁成分を混合させた。
【０１１８】
図３７は、心磁成分と特定部位の等価電流双極子との混合波形を説明する図である。等価電流双極子のIntensityは１０［ｎＡ・ｍ］に設定し、誘発磁場を作成した。心磁成分の振幅は、誘発磁場に比べて大きいので、０．１倍し、混合した。
【０１１９】
図３８は、１次聴覚野、１次視覚野に電流双極子があった時の等磁場線図を説明する図である。それぞれの電流双極子が作る等磁場線図を表している。
【０１２０】
作成した磁場データを用いて、等価電流双極子推定を行った。推定に使用したチャネルを図３９に示す。
【０１２１】
図３９は、解析に用いたチャネルを説明する図である。青枠は１次聴覚野近傍の等価電流双極子推定に用いた。赤枠は１次視覚野近傍の等価電流双極子推定に用いた。
【０１２２】
図４０、図４１に混合した磁場データに対し、等価電流双極子推定を行った結果を示す。
【０１２３】
図４０は、心磁成分が等価電流双極子に与える影響を説明する図である。上図において、色が付いた縦の棒はIntensityを表し、ドットはＧＯＦを表す。それぞれ、下のＭＲＩ断層写真上に示された等価電流双極子の位置と色が対応している。１次聴覚野近傍の電流双極子が心磁成分によって位置が下に移動している。Intensityにも影響を及ぼしている。
【０１２４】
図４１は、心磁成分が等価電流双極子に与える影響を説明する図である。上図において、色が付いた縦の棒はIntensityを表し、ドットはGOFを表す。それぞれ、下のＭＲＩ断層写真上に示された等価電流双極子の位置と色が対応している。１次視覚野近傍の等価電流双極子が心磁成分によって位置が大幅に移動している。また、Intensityも大きく変化している。図４０の１次聴覚野近傍より大きな影響を受けている。
【０１２５】
心磁成分が部位によって与える影響は異なっていた。図４０の等価電流双極子の位置は、主に上下に移動するのに対し、図４１は大きく脳幹や脳以外の場所に移動した。
【０１２６】
３．５．考察
本手法により、被験者４人から心磁成分の抽出を行うことが出来た。計測された心磁波形の違いは、椅子の座り方、被験者の心臓の位置、ＳＱＵＩＤからの心臓までの距離が異なる為だと考えられる。
【０１２７】
運動負荷は、心電のＱＲＳ時間、Ｑ−Ｔ間隔の短縮になって現れ、精神負荷は、Ｒ−Ｒ間隔の変動となって現れるという報告がある[8][9]。ＭＥＧ実験中は、運動を行う作業はしないので、心磁波形は相似になると考えていた。時間を分けて心磁波形を求めた結果（図３４）から、ほとんど心磁波形に変化は見られなかったので、本手法がＭＥＧ実験中でも使用できると考える。
【０１２８】
脳磁場解析を行う時、オフラインでフィルタ処理を行う。しかし、心磁の周波数成分を解析したところ、幅広く心磁成分は分布していた。心磁成分を除去する場合、通常用いられるフィルタ処理では、誘発磁場応答の成分も減衰させてしまい、除去するのは困難ということが分かった。その為に、心磁成分を除去する価値は十分にあると考えられる。
【０１２９】
抽出した心磁成分に対し、等価電流双極子を行った。この結果は、推定アルゴリズムによって異なると思うが、QRSに相当する心磁波形は、脳の中央位置（脳幹近傍）に等価電流双極子が推定された。
【０１３０】
心磁と混合した磁場データを人工的に作成して、等価電流双極子推定を行ったところ、大きく等価電流双極子位置がずれた。１次聴覚野近傍の等価電流双極子は、主に下方向に、１次視覚野近傍の等価電流双極子は脳幹〜脳幹下の位置に推定された。等価電流双極子は、等磁場線図の境目が出来る位置に推定される。心臓磁場成分が作る等磁場線図とまとめたものを図４２に示す。
【０１３１】
図４２は、心磁成分が１次聴覚野近傍、１次視覚野近傍に与える等磁場線図である。磁場の吹出と吸い込みの位置に等価電流双極子は推定される（黄矢印）。１次聴覚野が作る等磁場線図に対し、心磁の等磁場線図が与える影響は、吹出と吸い込みが影響を及ぼす。１視覚野が作る等磁場線図に対し、心磁の吹出と吸い込みの境目が影響を及ぼす。
【０１３２】
１次聴覚野近傍の等価電流双極子が作る等磁場線図と心磁成分が作る等磁場線図を比較すると、心磁成分は主に吸い込みと吹出となり影響を与える。１次視覚野の場合は、心磁成分の等磁場線図と境目が似てしまっている。その為に、１次視覚野の等価電流双極子は心磁の影響を強く受けてしまい、脳幹近傍に推定されたと考えることが出来る。
【０１３３】
これらの結果は、誘発磁場より心磁が強い条件である。加算平均回数が増えれば増えるほど、心磁の振幅は小さくなるので影響が小さくなると考えられる。その為に、実際に行っているMEG実験の脳磁場データとは異なる。しかし、心磁成分の影響が等価電流双極子の位置により、異なることを示すことが出来た。
【０１３４】
４．心磁成分を除去するプログラム
既存の解析ソフトウェアは、同時に記録した心電波形と加算平均により心磁成分の抽出を行うことができる。しかし、計測した脳磁場データから、心磁成分の除去は行えない。本研究で得られた知見と有効性を確認する為に、心磁成分を抽出、除去を行うプログラムを開発した（図４３）。プログラムの使用方法については付録で解説する。ここでは、プログラムの処理と心磁波形の除去結果について示す。
【０１３５】
図４３は、心磁成分を除去するプログラムのスクリーンショット（2008年1月現在）を示す図である。開発環境には、Microsoft Visual C++ 2005(Microsoft Corporation)を用いた。ＧＵＩ開発には、ＭＦＣ（Microsoft Foundation Class）を用いた。プログラム名は心磁を引くという意味をこめ「SubMCG」という名前にした。
【０１３６】
４．１．製作したプログラムのフローチャート
製作したプログラムのフローチャートを図４４に示す。本プログラムは、同時計測された心電波形から、閾値を越えた時刻をトリガに設定し加算平均波形を作り、心磁成分を抽出する。得られた、心磁成分から、トリガのタイミングにあわせて、ＳＱＵＩＤセンサから心磁成分を減算する。
【０１３７】
図４４は、プログラムのフローチャートである。プログラムは、大きく分けて４つの作業に分けられる。（１）心電波形の前処理、（２）心電波形からトリガ抽出、（３）心磁成分の抽出の為の加算平均、（４）ＳＱＵＩＤセンサの計測データから心磁成分を除去する。
【０１３８】
（１）心電波形の前処理
心電波形へフィルタ処理する。全ＳＱＵＩＤチャネルについてチャネル間平均値を求める。心電波形の最大値を求める。全サンプルの探索が終わるまで、チャネル間平均値を求め、心電波形の最大値を求める処理を繰り返す。
【０１３９】
（２）心電波形からトリガ抽出
心電波形の最大値から閾値を求める。心電波形と閾値から心磁成分抽出トリガを検出する。全サンプルの探索が終わるまで、心電波形と閾値から心磁成分抽出トリガを検出する処理を繰り返す。
【０１４０】
（３）心磁成分の抽出の為の加算平均
トリガ時間の前後に大きな変化があれば、加算平均のステップをスキップし、トリガ時間の前後に大きな変化がなければ、加算平均を行う。全サンプルの探索が終わるまでこれを繰り返す。
【０１４１】
（４）ＳＱＵＩＤセンサの計測データから心磁成分を減算
心磁成分抽出トリガから各ＳＱＵＩＤ計測データから心磁成分を減算する。全サンプルの探索が終わるまでこれを繰り返す。
【０１４２】
４．２．プログラムを用いた心磁成分の除去結果
製作したプログラムで、心磁成分を除去した結果を示す（図４５）。
【０１４３】
図４５は、本プログラムを用いた心磁成分の除去結果（プログラムのスクリーンショット）を説明する図である。心電波形に同期した心磁成分が計測されている。下のウィンドウには、心磁成分を除去した結果と、抽出した心磁成分を重ねて描いている。
【０１４４】
５．心磁成分除去システムの評価
心磁成分を除去するシステムを製作した。このシステムの評価を行うためにＭＥＧ実験を行った。
【０１４５】
５．１．実験条件
１次聴覚野、１次視覚野の誘発脳磁場応答が計測できる実験を設計した。
聴覚実験は、参考文献[4]の実験システムと音響刺激装置[19]を用いた。音刺激は、被験者の左耳に図４６のタイミングで呈示した。呈示した音刺激は、１０００［Ｈｚ］の正弦波で、音圧は６０［ｄＢ（Ａ）］である。音圧の測定には（株）小野測器のLA-5570（校正済）を用いた。計測はＡ特性（JIS C1505 Type1）を用いた。バックグラウンドノイズは１７〜２５ｄＢ（Ａ）であった。
【０１４６】
図４６は、呈示した音刺激を説明する図である。音刺激は２００［ｍｓ］で、立ち上がり立下り時間はそれぞれ１０［ｍｓ］であった。呈示間隔は、７００［ｍｓ］〜１０００［ｍｓ］で変化する。
【０１４７】
視覚実験は、参考文献[20]の実験システムを用いた。視覚刺激には、Ｃｏｓ関数で近似したガボールを被験者に提示した（図４７）。輝度測定には、株式会社ミノルタ（現コニカミノルタ）の分光放射輝度計CS-1000Aを用いて、輝度を測定した（図４７）。
【０１４８】
図４７は、呈示した視覚刺激を説明する図である。バックグラウンドの輝度は、１５．１［ｃｄ／ｍ２］、ガボール刺激の平均輝度は、２５．２［ｃｄ／ｍ２］であった。
【０１４９】
５．２．解析方法
脳磁場データ解析には、Meg Laboratory（横河電機、イーグルテクノロジー製）を使用した。被験者に提示した刺激をトリガにし、加算平均を行った。加算平均時、眼球運動や極端に大きな波形（２［ｐＴ］以上）は除いた。加算平均したデータに対し、カットオフ周波数が２−５０［Ｈｚ］のバンドパスフィルタをかけた。そのデータに対し等価電流双極子推定を行った。推定アルゴリズムはSarvas 法を用い、電流双極子モデルにはMoving Dipole（位置が動くことを想定したモデル）を用いた。結果で用いている電流双極子位置のＸ、Ｙ、Ｚ軸はそれぞれＭＲＩ座標系である（図４８）。
【０１５０】
図４８は、ＭＲＩ座標系を説明する図である。
【０１５１】
５．３．本手法の評価方法
心磁成分の除去効果を示す為に、心磁成分を除去した脳磁場データと、除去していない脳磁場データを用いて、等価電流双極子推定を行った（図４９）。また、加算平均回数を変えて比較した。加算平均回数は、２０回、５０回、１００回、３００回以上である。
【０１５２】
図４９は、心磁成分除去の評価方法を説明する図である。加算平均回数を同じに、心磁成分除去と未除去のデータの等価電流双極子位置の差を求めた。
【０１５３】
５．４．実験結果
３００回以上の刺激を呈示した脳磁場データに対して、等価電流双極子推定を行った結果を示す（図５０、図５１)。１次聴覚野、１次視覚野の推定結果は、位置、潜時、共に知られている知見と一致する。局在した時間を用いて、心磁成分の除去の効果を比較した。
【０１５４】
図５０は、聴覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。刺激開始から７０〜１３０［ｍｓ］の間で、1次聴覚野に等価電流双極子が局在した。
【０１５５】
図５１は、視覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。刺激開始から１１５〜１３５［ｍｓ］の間で、１次視覚野に等価電流双極子が局在した。
【０１５６】
心磁成分を除去した脳磁場データと比較を行った結果を示す。1次聴覚野での心磁成分除去の比較結果を図５２〜図５５に示す。
【０１５７】
図５２は、１次聴覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。Ｚ軸が大きく変化した。
【０１５８】
図５３は、１次聴覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。比較元のデータと差が少なくなってきた。Ｚ軸に差が出ている。
【０１５９】
図５４は、１次聴覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。Ｘ、Ｙ軸は図５３と同じく差がなくなってきた。Ｚ軸は時間共に変化していた。
【０１６０】
図５５は、１次聴覚野における等価電流双極子位置の差を説明する図である。加算平均は比較するデータと同じ回数行った。図５２〜図５４と比較して、Ｘ、Ｙ、Ｚ軸は平らになっていた。
【０１６１】
１次視覚野での心磁成分除去の比較結果を図５６〜図５９に示す。推定結果で４０［ｍｍ］以上のずれは省いた。
【０１６２】
図５６は、１次視覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。Ｘ、Ｙ、Ｚ軸が１０［ｍｍ］以上ずれていた。
【０１６３】
図５７は、１次視覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。図５６と比べて、時間と共に位置が変化しなくなってきた。Ｘ、Ｙ、Ｚ軸は１０〜２０［ｍｍ］間でずれている。
【０１６４】
図５８は、１次視覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。Ｚ軸が大きく変化した。Ｘ、Ｙ軸は元のデータに近づいた。
【０１６５】
図５９は、１次視覚野における心磁成分を除去したデータを用いた等価電流双極子位置の差を説明する図である。加算平均は比較するデータと同じ回数行った。Ｘ、Ｚ軸は比較元のデータと比べて差がなくなってきた。Ｙ軸は大きくずれていた。
【０１６６】
ＭＲＩ座標系のＸＹＺ軸の標準偏差の平均値を求め、加算平均回数の関係をまとめた（図６０)。
【０１６７】
図６０は、加算平均回数とＭＲＩ座標系ＸＹＺ軸の標準偏差を示す図である。加算平均回数が増えると、標準偏差の変化が一定になっていた。１次聴覚野では、加算平均回数５０回以降は大きな変化が無くなった。１次視覚野では加算平均回数１００回以降は大きな変化がなくなった。
【０１６８】
５．５．考察
本手法は、心磁成分のみを除去することが可能である。その為、加算平均は、ある程度必要である。１次視覚野の反応は、１次聴覚野の反応に比べて十分に小さい（図５０、図５１）ので、加算平均回数が少ない１次視覚野の結果は、２０〜３０［ｍｍ］ずれることになったと考える。
【０１６９】
１次聴覚野の結果は、Ｚ軸が大きく変化した。これは、３．４章で示唆した結果と同じになった。１次視覚野でもＺ、Ｙ軸に影響をおよぼす結果があった（図５８、図５９)。まとめると、心磁成分は図６１のように、推定位置に影響を与えることが分かった。
【０１７０】
図６１は、心磁成分の等磁場線図が等価電流双極子推定に与える影響を説明する図である。赤は磁場の吹出、緑は吸い込みをあらわす。青枠を使用して等価電流双極子推定を行うと、特にＺ軸方向に影響を特に受ける。赤枠を使用したチャネルは、吹出と吸い込みの境目の為に、心磁成分の影響を大きく受ける。
【０１７１】
視覚、聴覚の比較結果で加算平均回数が増すと、時間と共に差が変化しなくなり、一定の差だけが現れるようになった（図５５、図５９）。また、図６０から加算平均回数を増しても、一定の差だけが生じていることが分かる。この一定の差は、心磁成分だと考えられる。
【０１７２】
解析データは、２−５０［Ｈｚ］のバンドパスフィルタをかけているので、直流のような成分は減衰されている。この事から、心磁成分は基線をゆがませる一要因になっているということがいえる。通常用いられる基線補正は、計測データを直接変化させる為に、どこを基準に減算するかによって、大きく計測データをゆがませてしまう。また、基線が揺らいでいるデータに対しては有効ではない。本手法を用いることにより、基線をゆがませている一要因の心磁の除去と全区間に対しての基線補正をかけている事と同じになる。
【０１７３】
本手法を用いると、心磁が作る基線成分と推定位置への影響を除外できる。本実験では、１００回の脳磁場加算平均データを用いることにより、表２で表す程度の影響を除外できた。
（表２）
心磁成分を除去したことによる等価電流双極子の平均位置
実験の種類 MRI X軸[mm] MRI Y軸[mm] MRI Z軸[mm]
聴覚 0.1±1.6[mm] 2.4±2.4[mm] 1.5±3.2[mm]
視覚 4.6±2.8[mm] 2.0±2.0[mm] 16.4±5.6[mm]
【０１７４】
６．まとめ
本研究により、可能となった事を下記に示す。
シールドルーム内で実験のタスクを行っても、心拍を計測できる心電計の開発、
心電を用いて加算平均することにより、心磁成分を抽出するプログラムの開発、
計測されたSQUIDチャネルデータから心磁成分を除去するプログラムの開発。
【０１７５】
また、明らかになった事を下記に示す。
MEGでも心磁成分が計測され、フィルタ処理等では除去できない。
心磁成分が与える影響は、測定部位により異なることを示した。
頭部を上から見て、左右の位置はＭＲＩ座標系Ｚ軸に影響を及ぼし、真ん中の位置は、ＭＲＩ座標系Ｚ、Ｙ軸に影響を及ぼす。
心磁成分は加算平均を行うと、基線のゆらぎの一要因となってあらわれる。
本手法を用いることで、心磁成分を除去し、その効果示すことができた。
【０１７６】
７．参考文献
[1]樋口正法、磁束計支援システム , Med.Imag.Tech.Vol.1l No.2 June 1993.
[2]入戸野宏、堀忠雄、心理学研究における事象関連電位 (ERP) の利用、広島大学総合科学部紀要IV理系編、第26巻、pp15-31.
[3]小杉幸夫、武者利光（2000）、電子情報通信工学シリーズ 生体情報工学、森北出版株式会社、第1版、pp13-14.
[4]野村亮輔(2006)、Shepard Tone聴音時の脳磁場反応、平成18 年度金沢工業大学人間情報工学科卒業論文,pp7-25.
[5]小林哲生、栗城真也、多チャネル誘発脳磁界計測における自発脳磁界成分低減の検討、信学技報1996-06、pp7-12.
[6]池田思朗、村田昇、Independent Component Analysis を用いた MEG データの解析、電子情報通信学会技術研究報告、Vol.98、No.129、pp. 29-36.
[7]甲斐島 武、水野（松本） 由子、伊達 進、篠崎 和弘、下條 真司、MEGデータの独立成分分析（ICA）を用いた信号成分のクラスタリング解析、日本生磁気学会誌 特別号、Vol. 16、 No. 1、pp. 256-257.
[8]平野三千代、運動負荷後QRS時間およびQT間隔の回復時間と心電図異常との関係、岩手女子看護短期大学紀要、Vol.2(19950831)、pp15-21.
[9]山岸茂孝、津田良一、精神作業負荷心電図の周波数解析、第19回熊本県産学官技術交流会講演論文集、2005年1月、pp206-207
[10]トランジスタ技術編集部、上田智章(2006)、トランジスタ技術 2006年1月号、株式会社 CQ出版、43(1)、2006年1月、pp223-233.
[11]山口巖 監修、塚田啓二 編著、心磁図の読み方、株式会社コロナ社、2006年、第１版第１刷、pp10-12.
[12]山口武志、ベッドサイドのBasic Cardiology 心筋細胞の電気生理学―イオンチャネルから、心電図、不整脈へ―、メディカル・サイエンス・インターナショナル、２００６年、第１版第５刷、pp23-90.
[13]吉村正治(1970)、脈波判読の実際、株式会社中外医学社、第4版、pp1-31.
[14]DALE DUBIN,M.D.、図解心電図テキスト、株式会社 文光堂、2007年、第6版第３刷、p66.
[15]DALE DUBIN,M.D.、図解心電図テキスト、株式会社 文光堂、2007年、第6版第3刷、p36-54.
[16]イブ トーマス、中村 尚五、プラクティス デジタル信号処理、東京電機大学出版局、2005年、第１版４刷、pp132-154.
[17]村 上 正 浩、鈴 木 博 之、内 藤 茂 昭(2006)、日立心臓磁気計測システム MC-6400、Lab.Clin.Pract.、24(2)、pp132-139.
[18]トランジスタ技術編集部、山口晶大(2002)、トランジスタ技術 2002年10月号、株式会社 CQ出版,39(10)、2002年10月、pp173-186.
[19]冨田教幸(2005)、音階知覚に関連する脳磁場応答、平成17 年度金沢工業大学人間情報工学科卒業論文、pp41-56.
[20]山口和夫(2003)、Cognitive Timerに関連する脳磁場応答の計測、平成15年度金沢工業大学人間情報工学科卒業論文、pp40-50.
[21]山内光哉、心理・教育のための統計法<第２版>、株式会社サイエンス社、第２版 第８刷、pp68-78.
【０１７７】
８．付録
８．１．製作したプログラムについて
プログラムの使用方法について述べる。ＳＱＵＩＤチャネルの数、全チャネル数のデフォルト値は人間情報システム研究所のＭＥＧ計測装置にあわせている。
【０１７８】
図６２は、心磁成分除去プログラムの画面説明図である。
黄線の枠のOpen Binボタン、Searchボタン、SUB Preボタン、SUB Outボタンを順に押すことにより心磁成分が除去される。赤枠は、閾値の設定や読み込む心電チャネルを設定する。
【０１７９】
使用方法
１．ECG CHとフィルタを選ぶ。
２．ファイル読み込み(Open Binボタン)を行う。
３．Searchボタンを押すと心電波形ウインドウとトリガリストが更新される。
４．検出されたトリガを選択するか、波形の位置を選択する。すべて、正常にＥＣＧ波形のＲ波成分に赤線がかかっていた場合は５に進む。
トリガリストの数が極端に少ない場合、又は極端に多い場合は下記のことを試す。
Ｉ．心電波形が反転している場合は、反転をチェックしSearchボタンを押す。
ＩＩ．心電波形の閾値（赤線）が高すぎる、低すぎる場合はECG Coeffを設定する。本プログラムの閾値は、心電波形のPeak Threshold Ratioを求めている。ECG Coefの値を大きくすると閾値は高くなる。
ＩＩＩ．心電波形がＲ波ではない所に、トリガが設定されていた場合は、そのトリガを選択しDeleteボタンを押す。すると、そのトリガは削除される。
ＩＶ．読み込んだ波形が、心電ではなかった場合は、１からやり直す。
５．Sub Preボタンを押すと、心磁成分除去の為の前処理とMEG CHで選択したデータに対して、心磁成分を除去した結果を示せるようになる。
６．Sub Outボタンで心磁成分が除去されたファイルが出力される。
読み込んだファイルの名前が”abc.bin”だった場合は、”abc_MCG.con”というファイル名で、同じフォルダに出力されている。
【０１８０】
トリガリストのインポート、エクスポート
トリガリストのImport、ExportボタンはMEG Laboratoryのトリガのインポート、エクスポートに対応している。MEG Laboratoryを用いて、トリガをインポート、エクスポートすることも可能である。
【０１８１】
本プログラムで、読み込むことが出来る脳磁場データを図６３に示す。MEG Laboratoryは、図６３のような構造を持つバイナリファイルをエクスポート可能である。
【０１８２】
図６３は、プログラムが読み込むことの出来るバイナリファイル構造を説明する図である。人間情報システム研究所のＭＥＧ計測装置は、１サンプルあたり１９２ｃｈのデータを記録している。本プログラムもこの仕様にあわせた。１ｃｈ、１サンプルのデータはdouble型の８バイトである。
【０１８３】
本プログラムでは、計測した心電波形に対して、デジタルフィルタをかけて、基線の揺らぎとハムノイズを除去する。ＩＩＲフィルタでバンドパスフィルタを行うと、位相ひずみの為に、波形がひずんでしまう（図６４）。その為に、本プログラムでは、ＦＩＲ型のバンドパスフィルタとＩＩＲ型のバンドエリミネーションフィルタの２種類のフィルタを実装している。ＩＩＲフィルタのフィルタ係数の算出には、株式会社デジタルフィルターのDSPLinksを用いた。ＦＩＲフィルタのタップ係数の算出は、参考文献[14]を参考にして行った。
【０１８４】
図６４は、位相ひずみによる波形ゆがみを説明する図である。Raw Dataに対して、３−５０［Ｈｚ］のバンドパスフィルタをかける。ＩＩＲフィルタは、ゆがみ（図の※）が生じてしまっている。ＦＩＲフィルタは、ゆがみが生じない。
【０１８５】
本プログラムは、心電と同時記録することにより、心磁成分を除去することが出来る。その為に、心電を同時計測していないデータには適用できない。そこで、変化がはっきりとわかる、チャネルから、心電のトリガ抽出を試みた。波形の変化が顕著なチャネル（１０８ＣＨ）と、典型的な心磁波形からピアスンの積率相関係数[21]を求め、トリガ抽出を試みた。
【０１８６】
典型的な心磁波形は、被験者４人分の平均値を用いた（図６５）。変化の顕著なＱＲＳに相当する部分を用いた。
【０１８７】
図６５は、典型的な心磁波形を示す図である。計測した被験者４人分の心磁波形の平均値を求めて、変化の顕著な部分のみを用いた。
【０１８８】
図６６は、相関を用いてトリガの抽出を試みた結果を説明する図である。青線が心電トリガの時刻を表している。赤線は、１０８ｃｈと図６５の波形との相関を表す。ｒ＝０．８５以上の時をプロットした。
【０１８９】
図６６の結果から、心電より求めたトリガの時間と、相関より求めた時間が一致している部分もあるが、一致しない部分も出力された。単純に相関を求めるだけでは、トリガ抽出は出来なかった。何らかの拘束条件を用いて、探索する必要があると考える。例えば、ヒトの心電間隔（Ｒ−Ｒ間隔）は、安静時ほぼ一定になる特性を利用し除外するとよいかもしれない。他にも変化のはっきりとわかる他のチャネルを複数用いれば、ＳＱＵＩＤセンサから計測されるデータだけで、心電トリガの抽出が可能になるかもしれない。
【０１９０】
８．２．本手法を用いたハムノイズの除去
ハムノイズの原因は交流電源である。そこで、交流電源の周期を手に入れれば、本手法によりSQUIDセンサに混入するハムノイズのみ除去が出来るはずである。計測時に６０［Ｈｚ］のバンドエリミネーションフィルタをかけていない磁場データの周波数解析結果を示す（図６７）。
【０１９１】
図６７は、計測時に６０［Ｈｚ］のバンドエリミネーションフィルタをかけていない磁場データの周波数分布を示す図である。解析に使用したＳＱＵＩＤセンサ４ｃｈ分の平均した周波数分布を示す。６０［Ｈｚ］の成分とその高調波である１２０［Ｈｚ］、１８０［Ｈｚ］にピークが見られる(図の※)。
【０１９２】
オフラインにより、６０［Ｈｚ］のバンドエリミネーションフィルタをかけてハムノイズを除去した結果を図６８に、本手法を用いてハムノイズを除去した結果を図６９に示す。
【０１９３】
図６８は、バンドエリミネーションフィルタを用いたハムノイズ除去の結果を説明する図である。６０［Ｈｚ］成分のみが減衰している（赤矢印）。しかし、高調波成分は除去されない。
【０１９４】
図６９は、本手法を用いたハムノイズ除去の結果を説明する図である。６０［Ｈｚ］、高調波成分の１２０［Ｈｚ］、１８０［Ｈｚ］のピークが消えた。バンドエリミネーションフィルタを用いた結果と比較して、６０［Ｈｚ］の成分のみ減衰していない。
【０１９５】
ハムノイズの原因である交流電源と位相を完全に固定して除去している為に、６０［Ｈｚ］の脳磁場成分は減衰させることなく除去できる。その為に、ハムノイズ除去にも本手法が有効だと考えられる。
【０１９６】
なお、プログラムの開発にはMicrosoft CorporationのMicrosoft Visual C++ 2005を用いた。聴覚実験用の音刺激作成プログラムは、参考文献[19]のプログラムを使用した。
【０１９７】
以上、本発明を実施の形態をもとに説明した。実施の形態は例示であり、それらの各構成要素や各処理プロセスの組合せにいろいろな変形例が可能なこと、またそうした変形例も本発明の範囲にあることは当業者に理解されるところである。
【図面の簡単な説明】
【０１９８】
【図１】本実施の形態に係る脳磁場計測装置の構成図である。
【図２】ＭＥＧで計測される信号を説明する図である。
【図３】ＭＥＧで計測された生体現象によるノイズを説明する図である。
【図４】ＭＥＧで計測された環境が起因となっているノイズを説明する図である。
【図５】ランダムなノイズの自己相関を説明する図である。
【図６】時間相関のあるノイズを説明する図である。
【図７】MEG実験中の風景を説明する図である。
【図８】本研究が提案する心磁成分を除去する手法を説明する図である。
【図９】基線補正を説明する図である。
【図１０】心臓の各部位の名前と特徴を説明する図である。
【図１１】実験タスクの例を説明する図である。
【図１２】脈波観測装置を説明する図である。
【図１３】脈波観測装置の回路図である。
【図１４】製作した脈波観測装置を説明する図である。
【図１５】標準１２誘導の電極貼り付け位置を説明する図である。
【図１６】典型的な心電波形を説明する図である。
【図１７】製作した心電計の回路図である。
【図１８】製作した心電計を説明する図である。
【図１９】製作した電極リードと日本光電製の電極リードを説明する図である。
【図２０】空気動圧センサを用いた脈波計測の結果を説明する図である。
【図２１】光センサを用いた脈波計測の結果を説明する図である。
【図２２】フィルタをかける前の心電波形とバンドパスフィルタ後の心電波形を説明する図である。
【図２３】胸部誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。
【図２４】肢誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。
【図２５】心電計測中に電気刺激、ボタン押し動作をしたときの心電波形を説明する図である。
【図２６】光センサを用いた脈波観測装置と製作した心電計の心拍計測を説明する図である。
【図２７】電極の貼り付け位置を説明する図である。
【図２８】全頭型１６０チャンネル脳磁計を説明する図である。
【図２９】加算平均によって抽出した心磁波形と等磁場線図を説明する図である。
【図３０】各ＳＱＵＩＤセンサに含まれる心磁成分を説明する図である。
【図３１】心電波形と計測磁場データとの比較を示す図である。
【図３２】心電波形と計測した心磁図の比較を示す図である。
【図３３】被験者３人の心磁成分抽出結果を説明する図である。
【図３４】心磁波形の相似を説明する図である。
【図３５】心磁波形のパワースペクトル分布を説明する図である。
【図３６】抽出成分に対して等価電流双極子推定をおこなった結果を説明する図である。
【図３７】心磁成分と特定部位の等価電流双極子との混合波形を説明する図である。
【図３８】１次聴覚野、１次視覚野に電流双極子があった時の等磁場線図を説明する図である。
【図３９】解析に用いたチャネルを説明する図である。
【図４０】心磁成分が等価電流双極子に与える影響を説明する図である。
【図４１】心磁成分が等価電流双極子に与える影響を説明する図である。
【図４２】心磁成分が１次聴覚野近傍、１次視覚野近傍に与える等磁場線図である。
【図４３】心磁成分を除去するプログラムのスクリーンショットを示す図である。
【図４４】プログラムのフローチャートである。
【図４５】本プログラムを用いた心磁成分の除去結果（プログラムのスクリーンショット）を説明する図である。
【図４６】呈示した音刺激を説明する図である。
【図４７】呈示した視覚刺激を説明する図である。
【図４８】ＭＲＩ座標系を説明する図である。
【図４９】心磁成分除去の評価方法を説明する図である。
【図５０】聴覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。
【図５１】視覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。
【図５２】１次聴覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。
【図５３】１次聴覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。
【図５４】１次聴覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。
【図５５】１次聴覚野における等価電流双極子位置の差を説明する図である。
【図５６】１次視覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。
【図５７】１次視覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。
【図５８】１次視覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。
【図５９】１次視覚野における心磁成分を除去したデータを用いた等価電流双極子位置の差を説明する図である。
【図６０】加算平均回数とＭＲＩ座標系ＸＹＺ軸の標準偏差を示す図である。
【図６１】心磁成分の等磁場線図が等価電流双極子推定に与える影響を説明する図である。
【図６２】心磁成分除去プログラムの画面説明図である。
【図６３】プログラムが読み込むことの出来るバイナリファイル構造を説明する図である。
【図６４】位相ひずみによる波形ゆがみを説明する図である。
【図６５】典型的な心磁波形を示す図である。
【図６６】相関を用いてトリガの抽出を試みた結果を説明する図である。
【図６７】計測時に６０Ｈｚのバンドエリミネーションフィルタをかけていない磁場データの周波数分布を示す図である。
【図６８】バンドエリミネーションフィルタを用いたハムノイズ除去の結果を説明する図である。
【図６９】本手法を用いたハムノイズ除去の結果を説明する図である。
【符号の説明】
【０１９９】
１０ 生体信号計測部、 ２０ ノイズ信号抽出部、 ３０ ノイズ発生時刻取得部、 ４０ ノイズ除去処理部、 ５０ 生体計測信号記憶部、 ６０ ノイズ信号記憶部、 １００ 脳磁場計測装置。
【技術分野】
【０００１】
この発明は、生体計測データからノイズを除去する方法、ならびにその方法を利用可能な生体情報計測装置および脳磁場計測装置に関する。
【背景技術】
【０００２】
生体計測データには必ず、取得したい信号とノイズ成分が混合されている。たとえば、人間の脳磁場を測定すると、計測データには、脳磁場の信号以外に、生体現象による心電、商用電源から発生するハムノイズ、およびガウスノイズが混在している。これらのノイズのうち、心電とハムノイズは時間相関のあるノイズであり、ガウスノイズは時間相関のないノイズである。生体計測データから所望の信号を得るためにはノイズ除去が不可欠である。
【０００３】
従来行われてきたノイズ除去の方法として、計測を繰り返して加算平均を取ることで、ノイズ成分を除去する方法や、フィルタリング処理を施す方法、独立成分分析（Independent component analysis、ＩＣＲ）によって計測対象の信号と独立なノイズ成分を分離する方法などがある。
【発明の開示】
【発明が解決しようとする課題】
【０００４】
急峻な成分を持つ心電や心磁などの生体現象ノイズは、加算平均やフィルタ処理を行っても十分に減衰しない。また、ハムノイズの除去にはバンドエリミネーションフィルタ（ＢＥＦ）が用いられるが、ＢＥＦは、商用電源の周波数の成分を全て減衰させるため、元の信号成分も減衰してしまう。また、ＢＥＦでは、商用電源の高調波成分を減衰させることができず、ノイズとして残ってしまう。
【０００５】
心電や心磁などの生体現象ノイズは、測定対象の脳磁場の信号とは独立であるから、独立成分分析をすれば原理的には分離することが可能であが、独立成分分析は、勾配法などにより繰り返し探索をする必要があり、分析には非常に長い時間がかかる。また、独立成分分析の適用条件が厳しいため、実際の測定が条件をすべて満たすことは難しく、実用的ではない。
【０００６】
本発明はこうした課題に鑑みてなされたものであり、その目的は、生体計測データから効率良くノイズを除去する方法、ならびにその方法を利用可能な生体情報計測装置および脳磁場計測装置を提供することにある。
【課題を解決するための手段】
【０００７】
上記課題を解決するために、本発明のある態様のノイズ除去方法は、生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとを含む。
【０００８】
「時間相関性のあるノイズ信号の発生時刻に関する情報」は、ノイズ信号の発生時刻のタイムスタンプ情報であってもよく、ノイズ信号の発生時刻（トリガ）を求めるための基礎となるデータであってもよい。たとえば、ノイズ信号の大まかな波形データがあれば、ノイズ信号の値が所定の閾値を超える立ち上がり時点を検出したり、ノイズ波形のピークポイントを検出することで、ノイズ信号の発生時刻（トリガ）を求めることができる。
【０００９】
本発明の別の態様は、生体情報計測装置である。この装置は、生体から発生する信号を計測する計測部と、前記計測部による生体計測信号の計測の際に発生する時間相関性のあるノイズ信号の発生時刻に関する情報を取得する取得部と、前記計測部により計測された生体計測信号を前記時間相関性のあるノイズ信号の発生時刻に関する情報とともに格納する記憶部と、前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するノイズ信号抽出部と、前記時間相関性のあるノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得するノイズ除去処理部とを含む。
【００１０】
本発明のさらに別の態様は、脳磁場計測装置である。この装置は、生体から発生する脳磁場信号を計測する計測部と、前記脳磁場信号の計測時の心拍の発生時刻に関する情報を取得する取得部と、前記計測部により計測された脳磁場信号を前記心拍の発生時刻に関する情報とともに格納する記憶部と、前記脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出するノイズ信号抽出部と、心拍の発生時刻に合わせて、抽出された前記ノイズ信号を前記脳磁場信号から減算することにより、前記ノイズ信号が除去された脳磁場信号を取得するノイズ除去処理部とを含む。
【００１１】
なお、以上の構成要素の任意の組合せ、本発明の表現を方法、装置、システム、コンピュータプログラム、データ構造、記録媒体などの間で変換したものもまた、本発明の態様として有効である。
【発明の効果】
【００１２】
本発明によれば、生体計測データからノイズを効率的に除去することができ、被験者の負担を大幅に軽減することができる。
【発明を実施するための最良の形態】
【００１３】
本発明の実施の形態の概要を説明する。実施の形態では、生体計測データから時間相関のあるノイズ除去を行う方法を提供する。生体計測データから、時間相関のあるノイズを抽出して除去する。時間相関のあるノイズは、ガウスノイズとは違って、加算平均だけでは除去されないことが多い。そこで、実施の形態では、時間相関のあるノイズをノイズ信号自身に同期させ、いったんノイズ信号を計測データから抽出し、抽出されたノイズ信号を計測データから減算する。
【００１４】
図６に時間相関のあるノイズの自己相関が示されている。左図は周期性のあるハムノイズの自己相関であり、１６ｍｓの周期があることがわかる。右図は心電波形の自己相関であり、約７００ｍｓ付近に１回目のピークがあり、ピーク値は１．０より小さくなっている。また、２回目のピークは約１４００ｍｓにあるが、ピーク値は１回目よりもさらに低減している。これは、心電が、ハムノイズのように正確な周期的ノイズではなく、時間的に揺らいでいるためである。
【００１５】
生体計測データには、計測対象の信号にノイズ成分が混合されている。図２に一例を示すように、生体信号には、時間相関のあるノイズとして、生体現象による心電と商用電源から発生するハムノイズが、時間相関のないノイズとしてガウスノイズが混合されている。
【００１６】
心電や心磁など急峻な成分をもつ生体現象ノイズはフィルタ処理を施しても十分に減衰しない。また、通常、ハムノイズ除去にはバンドエリミネーションフィルタ（ＢＥＦ）が用いられるが、ＢＥＦは商用電源の周波数を減衰させることができるが、その高調波を減衰させることはできない。図６８に示すように、ＢＥＦを用いると電源周波数である６０Ｈｚの成分すべてが減衰してしまい、また、※印で示す高調波が残ってしまう。
【００１７】
一般に周波数フィルタリングをかけると、ノイズ成分だけではなく、原信号の周波数成分も捨ててしまうことになる。図３５に示されるとおり、心磁ノイズは０．５Ｈｚ〜４０Ｈｚまでの周波数成分を含み、これは、計測対象の脳磁場（数Ｈｚ〜数十Ｈｚ）とスペクトルがオーバーラップしている。このため、周波数フィルタリングの手法では、心磁成分のみを選択的に除去することは原理的に不可能である。直観的には、心磁を除去すると、脳磁場も除去することになるからである。
【００１８】
そこで、本手法では、計測時に、時間相関のあるノイズの時刻情報（以下「トリガ」という）と計測データとを同時に記録する。計測データには、計測信号とトリガが含まれている。計測信号をノイズのトリガの周期の幅で区切り、区切られた計測信号をトリガのタイミングで同期させて、加算平均を行うと、時間相関のあるノイズのみが抽出される（図８の左）。計測信号を心拍の時刻情報に同期させて加算平均を取ると、心磁成分だけが強調されて抽出される。ここで、心拍は時間的にゆらぎがあることに配慮して、加算平均時間の前後に余裕をもたせ、加算平均を取る時間の範囲を一部重複させている。加算平均を取る時間の範囲を重複させる幅は、心拍トリガが時間軸方向にゆらぐ量に応じて決める。このように加算平均を取る時間の範囲に余裕をもたせておくことで、心拍トリガのタイミングが多少前後にずれても、心拍トリガの周期で区切られた計測信号の間で心拍トリガが一致するように時間軸方向にずらして加算平均を取ることができる。
【００１９】
心電信号が急峻に変化するタイミングをトリガとして検出するが、トリガの検出方法としては、心電信号の値が閾値を超える立ち上がりポイントを検出すればよく、実施の形態ではこの方法が用いられる。トリガの別の検出方法は、心電信号のピークを検出することである。いずれにしても大まかな心電波形さえ得られれば、その心電波形から急峻な立ち上がりポイントをトリガとして検出することで心拍の時刻情報を取得することができる。
【００２０】
次に、加算平均によって抽出したノイズ成分を計測信号から除去する（図８の右）。具体的には、同時記録した心拍の時刻（心拍トリガ）に合わせて抽出されたノイズを元データから減算する。これにより計測された信号から心磁ノイズが除去される。このとき、加算平均前後に時間的余裕を持たせているので、心磁ノイズのような時間的にゆらぎを有するノイズにも対応できる。
【００２１】
また、本手法をハムノイズ除去に適用するとＢＥＦとは異なり、図６９に示すように、電源周波数の信号成分を減衰させることなく、高調波成分を含めてハムノイズだけを除去することができる。
【００２２】
このように、本手法は、ハムノイズ等の周期的なノイズの除去、生体現象による時間的にゆらぎのあるノイズ（たとえば、心拍ノイズ）の除去の両方に対応することができる。
【００２３】
本手法によって時間相関のあるノイズを生体計測データからあらかじめ除去することができ、残りのガウスノイズなどは加算平均によって除去すればよいから、通常のノイズ除去に必要な加算平均回数を半分から４分の１程度に減らすことができる。そのため、被験者の負担を大きく軽減することができる。たとえば、被験者に音を１秒に１回聴かせて脳磁を計測する実験を、従来は１０００回くらい繰り返して加算平均を取っていたため、実験には２０分近くかかっていた。これが少なくとも半分から４分の１程度の加算平均回数で済むため、実験は５〜１０分程度の時間に短縮される。
【００２４】
本手法は、脳波計測やＭＥＧ計測におけるハムノイズ除去や心磁ノイズ除去に応用できることを実施の形態で詳細に説明する。しかし、これは一実施例に過ぎず、本発明のノイズ除去方法は、脳波やＭＥＧ計測に限らず、一般に生体計測データから時間相関性のあるノイズを除去することを目的とするあらゆる技術分野に適用することができる。
【００２５】
図１は、本発明の実施の形態に係る脳磁場計測装置１００の構成図である。脳磁場計測装置１００は、後述のＳＱＵＩＤセンサを用いた脳磁図（ＭＥＧ）測定装置であるが、ここでは、脳磁場を測定するための一般的な構成については省略し、ノイズ除去に関する構成について説明する。
【００２６】
生体信号計測部１０は、生体から発生する脳磁場信号を計測して、生体計測信号記憶部５０に記録する。ノイズ発生時刻取得部３０は、同一の生体から測定される心電波形をもとに心拍の発生時刻に関する情報を取得して生体計測信号記憶部５０に記録する。生体信号計測部１０による脳磁場信号の記録と、ノイズ発生時刻取得部３０による心拍の発生時刻の記録は同時に行われる。
【００２７】
ノイズ信号抽出部２０は、生体計測信号記憶部５０に記録された脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出し、抽出されたノイズ信号をノイズ信号記憶部６０に記憶させる。
【００２８】
ノイズ除去処理部４０は、心拍の発生時刻に合わせて、ノイズ信号記憶部６０に記憶されたノイズ信号を生体計測信号記憶部５０に記録された脳磁場信号から減算することにより、ノイズ信号が除去された脳磁場信号を取得し、ノイズ信号除去後の脳磁場信号を生体計測信号記憶部５０に書き戻す。
【００２９】
計測された脳磁場信号には、心拍が原因となって発生する心磁成分がノイズとして混在する。本手法では、心磁信号を直接測定することなく、心電波形が得られる心拍の発生時刻情報だけを用いて、計測された脳磁場信号を心拍の時刻に同期させて加算平均を取ることにより、計測された脳磁場信号から心磁信号を抽出することができる。抽出された心磁信号を計測された脳磁場信号から差し引けば、心磁信号を除去して脳磁場信号だけを取り出すことができる。
【００３０】
計測対象の生体信号は、脳磁場信号に限られない。計測対象は、脳波信号であってもよく、その他の生体から発生する任意の信号が計測対象となりうる。また、計測対象の生体信号にノイズとして混在する信号は、心電信号に限らず、時間相関のある任意の信号がノイズとなりうる。たとえば、計測装置の交流電源から発生する周期性のあるハムノイズがノイズ信号として生体計測信号に混在する場合にも同様のノイズ除去手法を適用できる。すなわち、ハムノイズ信号の発生時刻情報を取得し、発生時刻に同期させて生体計測信号を加算平均することでハムノイズ信号を抽出することができる。抽出されたハムノイズ信号を生体計測信号から差し引けば、ハムノイズのない生体計測信号が得られる。なお、交流電源から発生するハムノイズの発生時刻は、交流がマイナス方向からプラス方向にゼロを切るタイミングを検出することで得られる。この発生時刻情報の記録を生体信号の計測と同時に行えばよい。
【００３１】
以下、本発明の実施の形態を発明者の「心磁成分が脳磁場計測データに及ぼす影響」に関する研究内容に即して詳細に説明する。脳磁場計測装置１００の各構成の動作の詳細は、下記の研究内容を参照して理解される。
【００３２】
本研究は、ＭＥＧ計測において生体現象ノイズの心磁成分に注目した。心磁成分は、カットオフ周波数３［Ｈｚ］程度のフィルタ処理だけでは減衰が困難である事と、注目しているＥＣＤの位置によって、その影響が異なることが分かった。そこで、加算平均を行い心磁成分の抽出を行い、これを除去するシステムの開発を行った。聴覚と視覚のＭＥＧ実験を行い、本システムの評価を行った結果、良好な性能を有することが示された。
【００３３】
以下の説明の章立てを示す。
１章では、ＭＥＧ実験中に計測された様々な時間相関のあるノイズの概観を解説し、心磁成分除去について問題を設定した。また、本手法で提案する心磁成分の除去手法について解説する。
【００３４】
２章では、ＭＥＧで使用できる心拍計測システムについて解説する。本研究で行う心磁成分の除去手法は、心拍をＭＥＧ計測中に同時に記録しなくてはならない。そこで、ＭＥＧ計測中に問題の起こらない心拍計測装置の開発と既存の計測装置との評価を行った。
【００３５】
３章では、ＭＥＧによる心磁成分の抽出を行った結果について解説する。抽出した心磁成分の特徴を明らかにし、従来行われているフィルタ処理等では減衰が困難だということを示した。また、人工的なデータを用いて、推定部位によって心磁の影響が異なることを示唆した。
【００３６】
４章では、既存の脳磁図解析ソフトウェアでは出来ない本手法による心磁除去を実現するプログラムについて解説する。
【００３７】
５章では、１次聴覚野、１次視覚野の脳磁場応答が計測できる実験を設計し本手法の評価を行った。その結果、心磁成分が加算平均により、基線の揺らぎの一要因となっていることを示し、本手法の有効性を示すことが出来た。
【００３８】
１．序論
ＭＥＧ計測はノイズとの戦いである[1]。ＭＥＧ（Magneto Encephalo Graphy：脳磁図）は、大脳皮質に生じる微弱な電流から生じる磁場を多数のＳＱＵＩＤ（Superconducting Quantum Interference Devices：超伝導量子干渉素子）センサを用いて計測する装置である。誘発脳磁場は大変微弱である為に、磁気シールドなどで外部の磁気ノイズを減衰させて計測を行う。しかし、この磁気減衰効果は、生体から発生する磁気ノイズに対しては働かない。その為に、計測される脳磁場データには、目的以外の信号が多数含まれる（図２）。通常は計測された脳磁場データに対して、加算平均や信号処理などの手法を用いて、Ｓ／Ｎ比（Signal/Noise Ratio）を向上させる。
【００３９】
計測されるノイズには、大きく分けて生体現象が原因のものと計測環境が原因のものの二つがある（図３、図４）。
【００４０】
図２は、ＭＥＧで計測される信号を説明する図である。計測される信号は、ノイズ成分、心磁成分（生体現象によるノイズ）、ハムノイズ（環境に起因するノイズ）等が、線形に混合された磁場信号が計測されている。
【００４１】
図３は、ＭＥＧで計測された生体現象によるノイズを説明する図である。１６０個のＳＱＵＩＤセンサから計測された磁場信号を重ね合わせた磁場波形と等磁場線図。等磁場線図の緑は磁場の吸い込み（Sink）、赤は磁場の吹出（Source）を示す。等磁場線図は黄線の時刻のものである。ａ）の眼球運動はＳＱＵＩＤセンサに近い為に大きな磁場変化がみられる。ｅ）の呼吸は、非常にゆっくりとした変化である。また、等磁場線図も他とは大きく違う。
【００４２】
図４は、ＭＥＧで計測された環境が起因となっているノイズを説明する図である。ＳＱＵＩＤセンサから計測された磁場信号とパワースペクトル分布を示す。６０［Ｈｚ］、１２０［Ｈｚ］、１８０［Ｈｚ］というハムノイズの６０［Ｈｚ］の倍数の周波数においてピークが見られる。このピークは、６０［Ｈｚ］の交流電源から発生していると考えられる。
【００４３】
ノイズには、時間に対して相関を持たないものと時間に対して相関を持つものがある（図５、図６）。誘発磁場以外の信号が、時間に対してランダムに変動すれば（時間と無相関な信号であれば）、Ｎ回の加算平均で１／Ｎの振幅に減衰する。これに対して、時間に対して相関のあるノイズ（時間相関のあるノイズ）信号はランダムに変動しないので、加算平均と同期してしまい、いつまでも相殺されずに残る可能性がある。また、生体現象ノイズ信号は通常脳磁場の１０倍程度の大きな信号である事が多いので、加算平均回数が１０回程度でも脳磁場と同程度の影響を結果に及ぼす。そのため、時間相関のあるノイズは、できるだけ加算平均から除くようにしなければいけない。一方時間相関があることを利用して、ノイズ成分を抽出できる可能性がある。
【００４４】
図５は、ランダムなノイズの自己相関を説明する図である。計測によって生じるノイズは、時間に対して相関を持たない。その為に、加算平均を行うと減衰する。
【００４５】
図６は、時間相関のあるノイズを説明する図である。左図は周期性のあるハムノイズの自己相関である。赤丸で示すようにピークが周期的に一定の間隔で生じている。右図は、心電の自己相関を表す。Ｒ−Ｒ間隔の部分を赤線で示す。心臓の拍動は一定の間隔で生じるが、時間的に揺らぎを生じる。
【００４６】
ＭＥＧ実験を行うとき、図７のように被験者の動きを抑制する。頭部が動くことによって生じる生体現象ノイズを抑制する。また、注視してもらうことで眼球運動が起こりづらくなるような実験タスクの工夫を行う。それでも、呼吸や心臓から、生体現象ノイズは生じてしまう。しかし、呼吸や眼球運動の生体現象のノイズは、通常計測される磁場データより数倍波形が大きく変化するので、加算平均時に含めないようにすることが可能である。また、電源から生じるハムノイズは計測時のバンドエリミネーションフィルタ（ＢＥＦ）、計測後のフィルタ処理により、大幅に減衰させることが可能である。しかし、心磁成分は、一定の間隔で生じ続け、除外することは困難である。
【００４７】
図７は、MEG実験中の風景を説明する図である。デュワーと被験者の頭部との間を、弾性がある素材で埋めている。これにより、頭部の横の動きは抑制される。また、野村が製作した顎乗せ台[4]を用いて、頭部の上下の動きを抑制している。顎乗せ台を使用する前のマーカーコイルずれは最大で９．７２［ｍｍ］であり、顎乗せ台を使用した後のマーカーコイルのずれは最大で０．４０［ｍｍ］になったという報告がある。
【００４８】
ＰＣＡやＩＣＡといった信号処理の手法による、ノイズ除去も提案されている[5][6]。独立成分分析（ＩＣＡ）は複数の観測信号から統計的な独立性を元に、源信号に分離する手法である。ノイズと誘発磁場は無相関である為に、互いに分離が可能である。しかし、ＭＥＧは多数のＳＱＵＩＤセンサを用いて計測するために、多次元のデータになり、解析には多くの時間がかかる[7]。
【００４９】
本研究の目的は、生体現象ノイズの心磁成分に注目し、心磁成分を除去するシステムの開発と等価電流双極子推定への結果に及ぼす影響を調べることで、このシステムの性能評価を行うことである。
【００５０】
本研究で行う心磁成分の除去手法は、心拍をトリガにし、加算平均を用いて心磁成分の抽出を行い、除去を行うという手法である。心磁成分は、一定のリズムで発生する為に、心拍を同時に記録すれば、心磁成分のみを除去するのに有効だと考えた。
【００５１】
１．１．加算平均による減算手法
本研究で提案する心磁成分の除去手法は、心臓の動きを、脳磁場計測と同時に記録する。これをトリガに用いて加算平均を行い、心臓の拍動によって生じる心磁成分を抽出し、オフラインで計測データから減算を行うという手法である（図８）。
【００５２】
一般に運動負荷を与えると心電波形のＱＲＳ時間、ＱＴ間隔が短縮する[8]。複写や計算課題などの負荷を与えると、Ｒ−Ｒ間隔が変動し精神的負荷を反映するといわれている[9]。本研究では、実験中ほぼ安静である為に抽出される心磁波形がほぼ相似になると考えた。Ｒ−Ｒ間隔の変動に対しては、１周期分の心磁波形を求め、心拍トリガから除去する事により対応が出来る。
【００５３】
図８は、本研究が提案する心磁成分を除去する手法を説明する図である。まず心拍の時刻情報を元に、ＳＱＵＩＤで計測された信号に対し、加算平均を行う。加算平均された信号は、心磁に由来する磁場成分であると考えられる。その後、加算平均に用いた心拍の時刻を用いて、抽出した心磁成分で除去を行う。
【００５４】
本手法を実現する為に、ＭＥＧ計測に影響を及ぼさない心拍計測装置と、本手法を行うソフトウェアが必要になる。
【００５５】
１．２．基線補正
加算平均を行ったデータに対し基線補正をかける。基線補正は、なんらかの持続的な神経活動や基線の揺らぎの除去を行う。基線補正が行う操作を図９に示す。
【００５６】
図９は、基線補正を説明する図である。指定した時間区間の各ＣＨの平均値を求めて、１サンプルごと元データを減算する。刺激トリガの前から１００［ｍｓ］、全区間で行われる。この操作は、直接計測データを変化させる。
【００５７】
直接計測データを変化させる為に、平均を求める区間によって、本来得たい信号を消してしまう可能性がある。また、どの部分を使えばよいという基準が明確に決められていない。
【００５８】
１．３．心電と心磁について
心臓は心筋を収縮させることによって、血液を全身に循環させるポンプの機能をもつ臓器である。心臓には、左右に隔壁と弁膜で仕切られた心室と心房が備わっている。それ以外にも、心房と心室を構成する心筋を収縮コントロールするための電気的指令を出したり、伝えたりする組織を供えている（図１０）。
【００５９】
図１０は、心臓の各部位の名前と特徴を説明する図である。洞結節は心筋の収縮の一定のリズムを作り出す部位である。左右房の心筋の脱分極による伝播は三尖弁と僧房弁でとまり、左右室への伝播は房室結節のみに伝わる。心室の刺激伝導系はヒス束から左右室の心筋まではプルキンエ繊維抹消である。この間は非常に早く伝わる。
【００６０】
心筋細胞は、静止時に陰性に保たれ、興奮によって脱分極する。この際、Ｎａ＋が細胞内に流入する。細胞はこの流入したＮａ＋を積極的に排出する。３個のＮａ＋が排出されると、２個のＫ＋が流入する。この時、外向きに電流が生じる。心筋が収縮する時Ｃａ２＋が必要になる。心筋の脱分極直後は、Ｎａ＋が多いのでＮａ＋を排出し、Ｃａ２＋を取り込む（外向き電流が生じる）。Ｃａ２＋が多くなると、逆にＣａ２＋を排出し、Ｎａ＋を取り込む（内向き電流が生じる）。心筋が収縮後、次の脱分極に備えて、再分極が起こる。この時、細胞内のＫ＋を排出し、元の静止膜電位に戻ろうとする[12]。
【００６１】
このように、心筋細胞はＮａ＋により脱分極、Ｃａ２＋により収縮、Ｋ＋により再分極がおこわなれる。心筋の興奮伝導（脱分極の伝播）を説明する。静止時、細胞内は陰性に保たれているが、隣り合う心筋細胞により脱分極が伝播する。この時、心筋の細胞内は陽性になる。興奮伝導の波面に電流双極子を生じる[11]。
【００６２】
体表面に電位差が生じる仕組みを説明する。心筋の興奮、収縮に伴う電流があると抵抗性を持つ周囲の臓器に帰還電流が生じる。その為に、任意の点で電位差が計測される（参考文献[10]の図２）。
【００６３】
心電波形と心筋細胞の収縮タイミングを説明する。計測される心電波形は、心筋の収縮タイミングを表している。参考文献[10]の図３に示されるように、心電図は心臓の電気的な状態を示し、Ｐ波、Ｑ波、Ｒ波、Ｓ波、Ｔ波、Ｕ波の順に現れる。Ｐ−Ｑの間に、左右房の収縮が起こる。ＱＲＳの間には、ヒス束からプルキンエ繊維抹消への興奮伝導が起こる。つまり、心室筋の脱分極のはじまりと終わりを表す。Ｓ−Ｔは、次の収縮の為に再分極している。
【００６４】
心磁の電流源は、興奮伝導に伴う電流双極子だと考えられる。心筋の興奮の伝播が心磁となって計測されると考えられる。
【００６５】
２．ＭＥＧで使用できる心拍計測システムについて
心臓は、心筋の収縮によって全身へ血液を運ぶ臓器である。この変化は、血管内の容積変化（脈波）と心筋の筋収縮によって生じる体表面の電位変化（心電図）によって計測できる。本研究の心磁除去の為に、心拍と同期した信号を、ＭＥＧ計測中に同時に記録しなくてはいけない。そこで、ＭＥＧ計測と同時に計測しても、ノイズや実験タスクに問題の発生しない心拍計測装置について検討を行った。
【００６６】
２．１．心拍を計測するシステムの評価
ＭＥＧ実験では、体性感覚用の刺激装置に電気刺激装置がある。また、実験タスクによっては、被験者にボタン押し動作を行う事がある。このような刺激や動作の時に、本件急で開発する心拍を計測するシステムが影響をうけないようにしなければいけない（図１１）。
【００６７】
図１１は、実験タスクの例を説明する図である。体性感覚用の刺激装置に電気刺激装を用いるものがある。実験タスクによっては、ボタン押し動作を行うことがある。これらの実験中に、心拍計測するシステムに影響を与えないようにしなければならない。
【００６８】
２．２．脈波を用いた心拍の計測
脈波とは、心拍動によって生ずる動脈系圧波動の伝達である[13]。この伝達を計測する方法に、空気動圧センサを用いて変動を計測する方法とヘモグロビンの赤外光における吸収量が違うことを利用する方法がある。
【００６９】
空気動圧センサを用いた脈波を計測する装置に（株）エム・アイ・ラボ製の脈波観測装置（図１２）を用いた。赤外線を用いた脈波観測装置は、赤外線ＬＥＤと光センサ（フォトダイオード）で構成した心拍計測装置を製作した（図１３は回路図、図１４は製作した装置）。
【００７０】
図１２は、脈波観測装置（空気圧の変化で捉える）を説明する図である。脈波観測装置は空気袋のわずかな空気圧の変動を捉え、電圧に変換し出力している。使用時は、空気袋を被験者の指にサージカルテープで密着させて使用する。このときチューブが空気を伝わる遅延と心臓から指までの血流の遅延を考慮しなくてはいけない。
【００７１】
図１３は、脈波観測装置（赤外線を使用した装置）の回路図である。赤外線LEDとフォトダイオードの間に指を挟むことにより、血流の変化を計測することができる。
【００７２】
図１４は、製作した脈波観測装置を説明する図である。赤外線は、人の目で見ることが出来ない。左図で、赤外ＬＥＤ部分が光っている。これは、デジタルカメラの光センサが赤外光に感受性をもつため記録される。右は、脈波の測定風景を示す。赤外線LEDとフォトダイオードを隠すように指を置いて計測を行う。
【００７３】
２．３．心電計の製作
心筋の細胞が興奮することによって生じる電場は、ほぼ、瞬間的に周辺の臓器を経由して体表面まで到達する。心電計は、この電場変化を捉える装置である。臨床目的の心電図は、標準１２誘導（肢誘導と胸部誘導）と呼ばれる電極の貼り付けを行う[15]（図１５）。本研究では、心電を臨床目的に使用するわけではないので、単純に心筋の収縮によって生じるＱＲＳの電場変化を捉えて増幅すればよい（図１６）。製作した心電計の回路図と装置概観を図１７、図１８に示す。
【００７４】
図１５は、標準１２誘導の電極貼り付け位置を説明する図である。肢誘導と胸部誘導の電極の取り付け位置を示す。この電極貼り付け位置は、臨床目的で使用されている。
【００７５】
図１６は、典型的な心電波形を説明する図である。心筋の細胞が、脱分極すると、心筋は収縮をする。Ｐ波は、心房の収縮を表す。ＱＲＳ波は、心室の脱分極の始まりから終わり、つまり心室の収縮を表し、Ｔ波で心筋の収縮を終えることを表している。詳細は１．３で解説している。
【００７６】
図１７は、製作した心電計の回路図である。ＡＤ６２０は計装アンプである。高い入力インピーダンスを持った差動増幅専用のオペアンプである。ＶＲで、ゲインを１〜１００００倍まで増幅できる。ＲＣで構成された０．５［Ｈｚ］のハイパスフィルタと４０［Ｈｚ］のローパスフィルタをかけた後、非反転増幅し出力している。
【００７７】
図１８は、製作した心電計を説明する図である。製作した心電計は、シールドルーム外に設置する。ケーブルをシールドルーム内まで伸ばし、被験者に電極を取り付ける。
【００７８】
被験者の体につける電極は、日本光電製ディスポ電極ビトロードF-150Sを用いた。電極は、実験前に磁石を用いて、磁性体ではないことを確認している。市販の製品は、ディスポ電極と心電計をつなぐ電極リードにステンレスが用いられており、ＭＥＧ実験に使用できないことが分かった。その為、磁性体を用いていない電極リードを自作した。導線は、銅を用い、ビトロードとの接続部分には、電子工作用に用いられるゴム足を加工して用いた（図１９）。
【００７９】
図１９は、製作した電極リード（左）と日本光電製の電極リード（右）を説明する図である。自作した電極リードは銅線部分がディスポ電極に接触するように出来ている。日本光電製の電極リード（使い捨て電極リード線BC−IXX：XXは色によって型番が異なる）はディスポ電極をバネの構造で密着させている。
【００８０】
計測した心電データに対して、デジタルフィルタをかける。このフィルタ処理は、計測中の被験者の動きによる基線の揺らぎ（〜２［Ｈｚ］）や、ハムノイズ（６０［Ｈｚ］の成分）を除去する為のものである。
【００８１】
２．４．脈波観測装置と製作した心電計の心拍計測結果
図１２と図１４の脈波観測装置を用いて、心拍を計測した結果を示す（図２０、図２１）。この計測の被験者は、共に、２秒間に３回ほど、振幅のピークが見られた。安静時の心拍の間隔が、正常な間隔が１秒間に０．８〜１．７回[14]なので、正常に心拍が計測されているものと考えられる。
【００８２】
図２０は、空気動圧センサを用いた脈波計測の結果を説明する図である。指を動かすと、大きく基線がゆらいだ。また、空気袋を触ると大きな変化になってしまい、脈波が正常に観測されない。
【００８３】
図２１は、光センサを用いた脈波計測の結果を説明する図である。基本的に指を光センサと赤外線LEDとの間に置き続けなければ計測が出来ない。その為に、ＭＥＧ実験中に、ボタンを押すという動作は出来ない。
【００８４】
計測した心電波形に、３−５０［Ｈｚ］のＦＩＲ型のバンドパスフィルタをかけた。これにより、被験者が動くことによって生じる心電波形の基線の揺らぎとハムノイズを除去できる（図２２）。この時、位相ひずみが起きないようにフィルタは直線位相特性を持ったものを使用した[16]。
【００８５】
図２２は、フィルタをかける前の心電波形とバンドパスフィルタ後の心電波形を説明する図である。細かく一定の周期で現れているのがハムノイズである。交流電源の周波数が６０［Ｈｚ］なので、その帯域を減衰するように３−５０［Ｈｚ］のバンドパスフィルタをかけた。
【００８６】
製作した心電計は２点間の電位差を求めて増幅することが出来る（予備実験の計測は、ノイズ低減の為に不感電極を心電計のＧＮＤに接続して計測を行っている。）。よって、どの部位に電極を貼り付けるかによって、波形は異なる。そこで、電極を頂点とみなし、三角形に心臓が収まるように電極をつけて計測を行った（図２３、図２４）。
【００８７】
図２３は、胸部誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。赤がプラス入力、青がマイナス入力、黒が不感電極（心電計ＧＮＤ）に接続している。
【００８８】
図２４は、肢誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。赤がプラス入力、青がマイナス入力、黒が不感電極（心電計ＧＮＤ）に接続している。全体的に、胸部誘導より、振幅が低くなった。これは、心臓からの距離が離れた為に、臓器などが抵抗成分となった為だと考えられる。
【００８９】
ＭＥＧ実験で使用する電気刺激装置とボタン押し動作をした時の結果を図２５に示す。心電波形に変化は見られなかった。
【００９０】
図２５は、心電計測中に電気刺激、ボタン押し動作をしたときの心電波形を説明する図である。赤い線は、それぞれの刺激の始まりのタイミングを表している。電気刺激は０．５［ｍＡ］を手首に出力している。ボタン押し動作は、人差し指でボタンを押してもらった。共に心電波形に変化は見られなかった。
【００９１】
製作した脈波観測装置（図１４）と心電計を用いて同時に心拍を計測した結果を図２６に示す。心電のＲ波と脈波のピーク値はおよそ２００［ｍｓ］の遅延が生じていた。この遅延は、被験者によって心臓から脈波の測定部位までの距離が異なる為に、変化する。
【００９２】
図２６は、光センサを用いた脈波観測装置と製作した心電計の心拍計測を説明する図である。脈波は心電に比べて遅延が生じていた。これは、血流の変化を捉えている為である。
【００９３】
２．５．考察
本研究では、MEG計測と心拍を同時に記録し、加算平均することにより心磁成分のみを抽出する。その為に、実験タスク中でも安定して心拍を計測する装置が必要になる。
脈波観測装置は、空気動圧センサを用いるもの、赤外線を用いるもの、それぞれに問題がある。実験タスクによって、ボタンを押すという動作をすることがある。図２０から、空気動圧センサを用いた脈波観測装置は、使用できない。光センサを使用する脈波観測装置は、センサをシールドルーム内に持ち込むことは出来ない。光センサの電流の変化が、ＳＱＵＩＤセンサに混入してしまうからである。また、光ファイバを用いることも考えたが、赤外線の波長の減衰が大きく、脈波成分が計測出来なかった。
【００９４】
製作した心電計は、電気刺激、ボタン押し動作を行っても、安定して心拍が計測できた。また、脈波とは違い、電場の変化を計測しているので、遅延が生じない。そこで、本研究では、心拍を計測する装置に心電計を用いて、心磁成分の抽出を行うことにした。
【００９５】
３．ＭＥＧによる心磁成分の抽出
生体から発生する磁気を用いて機能を計測する装置は多数存在する。脳磁図は文字通り、脳機能をターゲットに作られた計測装置である。それ以外にも、心臓の機能を調べる心磁図用の計測装置がある。本研究においては、ＭＥＧ計測装置を用いた心磁図の報告は無いので、心磁がどのようにＳＱＵＩＤセンサに計測されるか調べた。
【００９６】
参考文献[17]に示される心磁図計測装置を用いる。心臓の真上に、ＳＱＵＩＤセンサアレイを配置し、計測を行う。ＭＥＧ計測装置は、頭部を覆うようにＳＱＵＩＤセンサを配置している。
【００９７】
心磁図測定装置を用いて計測された心磁成分は、参考文献[17]に示されるように、ＳＱＵＩＤセンサ８ｃｈ分を重ね合わせた結果である。
【００９８】
３．１．実験方法
被験者に電極を図２７のように貼り付けた。また、心電計の出力をＭＥＧ計測装置のトリガチャンネルに接続し、心電を同時に記録した。被験者は２２−２３歳男性４名で行った。
【００９９】
図２７は、電極の貼り付け位置を説明する図である。心臓が電極を頂点に、三角形の中に捉えられるように、反転入力、非反転入力、不感電極を貼り付けた。
【０１００】
心磁成分の抽出には、同時に記録した心電図のＲ波部分を、トリガに加算平均を行った。解析に用いたソフトは、脳磁図用の解析ソフトウェアMEG Laboratory（横河電機、イーグルテクノロジー製）である。加算平均は２００回以上行った。通常の脳磁場の解析に用いられるフィルタ処理、基線補正処理は行わなかった。
【０１０１】
３．２．ＭＥＧ計測装置について
脳磁場計測と心磁成分の計測には全頭型１６０チャンネル脳磁計（イーグルテクノロジー製）を使用した。ＳＱＵＩＤセンサはGradio-meterでBaselineは５０．０［ｍｍ］であった。ノイズとシグナルのピーク値を用いて計算したＳ／Ｎ比は０．０３［ｄＢ］であった。被験者は磁気遮断室（ＭＳＲ：Magnetically Shielded Room）内に入ってもらう。脳磁場計測には、顎乗せ台を用い、デュワーと頭部の間を図７のように埋めて計測を行った。
【０１０２】
図２８は、全頭型１６０チャンネル脳磁計を説明する図である。被験者は椅子に座り、デュワー内に頭部を入れてもらい脳磁場計測を行う。
【０１０３】
全ての計測時に共通なパラメータを表１に示す。
（表１）
ＭＥＧ計測のパラメータ
サンプリング周波数［Ｈｚ］ １０００
Low Pass Filterのカットオフ周波数［Ｈｚ］ ２００
Band Elimination Filterの中心周波数［Ｈｚ］ ６０
【０１０４】
３．３．ＭＥＧによる心磁成分の抽出結果
心電波形をトリガに、加算平均を行った結果を示す。図２９は加算平均によって得られた波形と、図３０はSQUIDチャネルの波形変化を示す。
【０１０５】
図２９は、加算平均によって抽出した心磁波形と等磁場線図を説明する図である。等磁場線図の緑は磁場の吸い込み（Sink）、赤は磁場の吹出（Source）を示す。変化の激しい時刻（赤線、青線）の等磁場線図を示す。
【０１０６】
図３０は、各ＳＱＵＩＤセンサに含まれる心磁成分を説明する図である。各ＳＱＵＩＤセンサを大きく左右に分けると、左側で最も大きく変化するチャネルは１０８ｃｈ、右側で最も大きく変化するチャネルは７１ｃｈであった。
【０１０７】
心電をトリガに、図２９のような波形が得られたことから、心拍による誘発磁場だと考えられる。同時に記録した心電波形と図３０から分かった変化が著しいチャネルの計測磁場データとの比較を行ったものを図３１に示す。
【０１０８】
図３１は、心電波形と計測磁場データとの比較を示す図である。心電波形のＱＲＳ波にかけて、７１ｃｈ、１０８ｃｈに同期した変化が見られる。７１ｃｈから１０８ｃｈを減算するとより分かりやすい。
【０１０９】
図３２は、心電波形と計測した心磁図の比較を示す図である。Ｐ波、ＱＲＳ波、Ｔ波に対応した心磁成分が計測されている。ＱＲＳ波に相当する部分が大きく計測された。
【０１１０】
被験者３名の心磁成分抽出結果を図３３に示す。全被験者でＰ、ＱＲＳ、Ｔに相当する心磁波形が確認できた。
【０１１１】
図３３は、被験者３人の心磁成分抽出結果を説明する図である。等磁場線図は図２９の青線と同じ心磁成分の位置を示す。被験者によって、ＱＲＳ波に相当する心磁成分の波形が異なっていた。また、振幅のピーク値も被験者によって異なった。被験者４人計測して、振幅の絶対値を取った所、最大で８２０［ｆＴ］、最小で３３３［ｆＴ］、平均して４８５±２２９［ｆＴ］の振幅が計測された。
【０１１２】
心磁波形が時間的に変化するか調べた（図３４）。心磁成分の周波数成分を図３５に示す。
【０１１３】
図３４は、心磁波形の相似を説明する図である。図のＡ、Ｂ、Ｃは、それぞれ時間の区間を表す。約２００［ｓ］で３００回の加算平均を行った。心磁波形にほとんど変化は無かった。
【０１１４】
図３５は、心磁波形のパワースペクトル分布（１０８ｃｈのみ）を説明する図である。ＱＲＳの心磁成分が急峻な変化である為に、広い周波数範囲に成分が含まれている。図によるとほぼ１〜５０［Ｈｚ］までの成分が含まれている。１〜５０［Ｈｚ］の成分は、脳磁場解析にも含まれる帯域である。
【０１１５】
３．４．心磁成分が脳磁場計測データに与える影響
抽出した心磁成分に対し、等価電流双極子推定をおこなった結果を図３６に示す。
【０１１６】
図３６は、抽出成分に対して等価電流双極子推定をおこなった結果を説明する図である。推定位置の色とＧＯＦ、Intensityがそれぞれ対応している。心磁波形の振幅が最も高い時刻に、４２０［ｎＡ・ｍ］ものIntensityが計測された。推定位置は、脳中央に推定された。推定アルゴリズムはSarvas法を用いた。
【０１１７】
心磁成分が、誘発磁場応答に、どのように影響を与えているか調べる為に図３７のような磁場データを人工的に作成した。作成にはMEG Laboratoryを用いた。手動で右の１次聴覚野近傍と１次視覚野近傍に等価電流双極子を置いた（図３８）。それぞれの等価電流双極子が作る磁場データに対し、心磁成分を混合させた。
【０１１８】
図３７は、心磁成分と特定部位の等価電流双極子との混合波形を説明する図である。等価電流双極子のIntensityは１０［ｎＡ・ｍ］に設定し、誘発磁場を作成した。心磁成分の振幅は、誘発磁場に比べて大きいので、０．１倍し、混合した。
【０１１９】
図３８は、１次聴覚野、１次視覚野に電流双極子があった時の等磁場線図を説明する図である。それぞれの電流双極子が作る等磁場線図を表している。
【０１２０】
作成した磁場データを用いて、等価電流双極子推定を行った。推定に使用したチャネルを図３９に示す。
【０１２１】
図３９は、解析に用いたチャネルを説明する図である。青枠は１次聴覚野近傍の等価電流双極子推定に用いた。赤枠は１次視覚野近傍の等価電流双極子推定に用いた。
【０１２２】
図４０、図４１に混合した磁場データに対し、等価電流双極子推定を行った結果を示す。
【０１２３】
図４０は、心磁成分が等価電流双極子に与える影響を説明する図である。上図において、色が付いた縦の棒はIntensityを表し、ドットはＧＯＦを表す。それぞれ、下のＭＲＩ断層写真上に示された等価電流双極子の位置と色が対応している。１次聴覚野近傍の電流双極子が心磁成分によって位置が下に移動している。Intensityにも影響を及ぼしている。
【０１２４】
図４１は、心磁成分が等価電流双極子に与える影響を説明する図である。上図において、色が付いた縦の棒はIntensityを表し、ドットはGOFを表す。それぞれ、下のＭＲＩ断層写真上に示された等価電流双極子の位置と色が対応している。１次視覚野近傍の等価電流双極子が心磁成分によって位置が大幅に移動している。また、Intensityも大きく変化している。図４０の１次聴覚野近傍より大きな影響を受けている。
【０１２５】
心磁成分が部位によって与える影響は異なっていた。図４０の等価電流双極子の位置は、主に上下に移動するのに対し、図４１は大きく脳幹や脳以外の場所に移動した。
【０１２６】
３．５．考察
本手法により、被験者４人から心磁成分の抽出を行うことが出来た。計測された心磁波形の違いは、椅子の座り方、被験者の心臓の位置、ＳＱＵＩＤからの心臓までの距離が異なる為だと考えられる。
【０１２７】
運動負荷は、心電のＱＲＳ時間、Ｑ−Ｔ間隔の短縮になって現れ、精神負荷は、Ｒ−Ｒ間隔の変動となって現れるという報告がある[8][9]。ＭＥＧ実験中は、運動を行う作業はしないので、心磁波形は相似になると考えていた。時間を分けて心磁波形を求めた結果（図３４）から、ほとんど心磁波形に変化は見られなかったので、本手法がＭＥＧ実験中でも使用できると考える。
【０１２８】
脳磁場解析を行う時、オフラインでフィルタ処理を行う。しかし、心磁の周波数成分を解析したところ、幅広く心磁成分は分布していた。心磁成分を除去する場合、通常用いられるフィルタ処理では、誘発磁場応答の成分も減衰させてしまい、除去するのは困難ということが分かった。その為に、心磁成分を除去する価値は十分にあると考えられる。
【０１２９】
抽出した心磁成分に対し、等価電流双極子を行った。この結果は、推定アルゴリズムによって異なると思うが、QRSに相当する心磁波形は、脳の中央位置（脳幹近傍）に等価電流双極子が推定された。
【０１３０】
心磁と混合した磁場データを人工的に作成して、等価電流双極子推定を行ったところ、大きく等価電流双極子位置がずれた。１次聴覚野近傍の等価電流双極子は、主に下方向に、１次視覚野近傍の等価電流双極子は脳幹〜脳幹下の位置に推定された。等価電流双極子は、等磁場線図の境目が出来る位置に推定される。心臓磁場成分が作る等磁場線図とまとめたものを図４２に示す。
【０１３１】
図４２は、心磁成分が１次聴覚野近傍、１次視覚野近傍に与える等磁場線図である。磁場の吹出と吸い込みの位置に等価電流双極子は推定される（黄矢印）。１次聴覚野が作る等磁場線図に対し、心磁の等磁場線図が与える影響は、吹出と吸い込みが影響を及ぼす。１視覚野が作る等磁場線図に対し、心磁の吹出と吸い込みの境目が影響を及ぼす。
【０１３２】
１次聴覚野近傍の等価電流双極子が作る等磁場線図と心磁成分が作る等磁場線図を比較すると、心磁成分は主に吸い込みと吹出となり影響を与える。１次視覚野の場合は、心磁成分の等磁場線図と境目が似てしまっている。その為に、１次視覚野の等価電流双極子は心磁の影響を強く受けてしまい、脳幹近傍に推定されたと考えることが出来る。
【０１３３】
これらの結果は、誘発磁場より心磁が強い条件である。加算平均回数が増えれば増えるほど、心磁の振幅は小さくなるので影響が小さくなると考えられる。その為に、実際に行っているMEG実験の脳磁場データとは異なる。しかし、心磁成分の影響が等価電流双極子の位置により、異なることを示すことが出来た。
【０１３４】
４．心磁成分を除去するプログラム
既存の解析ソフトウェアは、同時に記録した心電波形と加算平均により心磁成分の抽出を行うことができる。しかし、計測した脳磁場データから、心磁成分の除去は行えない。本研究で得られた知見と有効性を確認する為に、心磁成分を抽出、除去を行うプログラムを開発した（図４３）。プログラムの使用方法については付録で解説する。ここでは、プログラムの処理と心磁波形の除去結果について示す。
【０１３５】
図４３は、心磁成分を除去するプログラムのスクリーンショット（2008年1月現在）を示す図である。開発環境には、Microsoft Visual C++ 2005(Microsoft Corporation)を用いた。ＧＵＩ開発には、ＭＦＣ（Microsoft Foundation Class）を用いた。プログラム名は心磁を引くという意味をこめ「SubMCG」という名前にした。
【０１３６】
４．１．製作したプログラムのフローチャート
製作したプログラムのフローチャートを図４４に示す。本プログラムは、同時計測された心電波形から、閾値を越えた時刻をトリガに設定し加算平均波形を作り、心磁成分を抽出する。得られた、心磁成分から、トリガのタイミングにあわせて、ＳＱＵＩＤセンサから心磁成分を減算する。
【０１３７】
図４４は、プログラムのフローチャートである。プログラムは、大きく分けて４つの作業に分けられる。（１）心電波形の前処理、（２）心電波形からトリガ抽出、（３）心磁成分の抽出の為の加算平均、（４）ＳＱＵＩＤセンサの計測データから心磁成分を除去する。
【０１３８】
（１）心電波形の前処理
心電波形へフィルタ処理する。全ＳＱＵＩＤチャネルについてチャネル間平均値を求める。心電波形の最大値を求める。全サンプルの探索が終わるまで、チャネル間平均値を求め、心電波形の最大値を求める処理を繰り返す。
【０１３９】
（２）心電波形からトリガ抽出
心電波形の最大値から閾値を求める。心電波形と閾値から心磁成分抽出トリガを検出する。全サンプルの探索が終わるまで、心電波形と閾値から心磁成分抽出トリガを検出する処理を繰り返す。
【０１４０】
（３）心磁成分の抽出の為の加算平均
トリガ時間の前後に大きな変化があれば、加算平均のステップをスキップし、トリガ時間の前後に大きな変化がなければ、加算平均を行う。全サンプルの探索が終わるまでこれを繰り返す。
【０１４１】
（４）ＳＱＵＩＤセンサの計測データから心磁成分を減算
心磁成分抽出トリガから各ＳＱＵＩＤ計測データから心磁成分を減算する。全サンプルの探索が終わるまでこれを繰り返す。
【０１４２】
４．２．プログラムを用いた心磁成分の除去結果
製作したプログラムで、心磁成分を除去した結果を示す（図４５）。
【０１４３】
図４５は、本プログラムを用いた心磁成分の除去結果（プログラムのスクリーンショット）を説明する図である。心電波形に同期した心磁成分が計測されている。下のウィンドウには、心磁成分を除去した結果と、抽出した心磁成分を重ねて描いている。
【０１４４】
５．心磁成分除去システムの評価
心磁成分を除去するシステムを製作した。このシステムの評価を行うためにＭＥＧ実験を行った。
【０１４５】
５．１．実験条件
１次聴覚野、１次視覚野の誘発脳磁場応答が計測できる実験を設計した。
聴覚実験は、参考文献[4]の実験システムと音響刺激装置[19]を用いた。音刺激は、被験者の左耳に図４６のタイミングで呈示した。呈示した音刺激は、１０００［Ｈｚ］の正弦波で、音圧は６０［ｄＢ（Ａ）］である。音圧の測定には（株）小野測器のLA-5570（校正済）を用いた。計測はＡ特性（JIS C1505 Type1）を用いた。バックグラウンドノイズは１７〜２５ｄＢ（Ａ）であった。
【０１４６】
図４６は、呈示した音刺激を説明する図である。音刺激は２００［ｍｓ］で、立ち上がり立下り時間はそれぞれ１０［ｍｓ］であった。呈示間隔は、７００［ｍｓ］〜１０００［ｍｓ］で変化する。
【０１４７】
視覚実験は、参考文献[20]の実験システムを用いた。視覚刺激には、Ｃｏｓ関数で近似したガボールを被験者に提示した（図４７）。輝度測定には、株式会社ミノルタ（現コニカミノルタ）の分光放射輝度計CS-1000Aを用いて、輝度を測定した（図４７）。
【０１４８】
図４７は、呈示した視覚刺激を説明する図である。バックグラウンドの輝度は、１５．１［ｃｄ／ｍ２］、ガボール刺激の平均輝度は、２５．２［ｃｄ／ｍ２］であった。
【０１４９】
５．２．解析方法
脳磁場データ解析には、Meg Laboratory（横河電機、イーグルテクノロジー製）を使用した。被験者に提示した刺激をトリガにし、加算平均を行った。加算平均時、眼球運動や極端に大きな波形（２［ｐＴ］以上）は除いた。加算平均したデータに対し、カットオフ周波数が２−５０［Ｈｚ］のバンドパスフィルタをかけた。そのデータに対し等価電流双極子推定を行った。推定アルゴリズムはSarvas 法を用い、電流双極子モデルにはMoving Dipole（位置が動くことを想定したモデル）を用いた。結果で用いている電流双極子位置のＸ、Ｙ、Ｚ軸はそれぞれＭＲＩ座標系である（図４８）。
【０１５０】
図４８は、ＭＲＩ座標系を説明する図である。
【０１５１】
５．３．本手法の評価方法
心磁成分の除去効果を示す為に、心磁成分を除去した脳磁場データと、除去していない脳磁場データを用いて、等価電流双極子推定を行った（図４９）。また、加算平均回数を変えて比較した。加算平均回数は、２０回、５０回、１００回、３００回以上である。
【０１５２】
図４９は、心磁成分除去の評価方法を説明する図である。加算平均回数を同じに、心磁成分除去と未除去のデータの等価電流双極子位置の差を求めた。
【０１５３】
５．４．実験結果
３００回以上の刺激を呈示した脳磁場データに対して、等価電流双極子推定を行った結果を示す（図５０、図５１)。１次聴覚野、１次視覚野の推定結果は、位置、潜時、共に知られている知見と一致する。局在した時間を用いて、心磁成分の除去の効果を比較した。
【０１５４】
図５０は、聴覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。刺激開始から７０〜１３０［ｍｓ］の間で、1次聴覚野に等価電流双極子が局在した。
【０１５５】
図５１は、視覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。刺激開始から１１５〜１３５［ｍｓ］の間で、１次視覚野に等価電流双極子が局在した。
【０１５６】
心磁成分を除去した脳磁場データと比較を行った結果を示す。1次聴覚野での心磁成分除去の比較結果を図５２〜図５５に示す。
【０１５７】
図５２は、１次聴覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。Ｚ軸が大きく変化した。
【０１５８】
図５３は、１次聴覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。比較元のデータと差が少なくなってきた。Ｚ軸に差が出ている。
【０１５９】
図５４は、１次聴覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。Ｘ、Ｙ軸は図５３と同じく差がなくなってきた。Ｚ軸は時間共に変化していた。
【０１６０】
図５５は、１次聴覚野における等価電流双極子位置の差を説明する図である。加算平均は比較するデータと同じ回数行った。図５２〜図５４と比較して、Ｘ、Ｙ、Ｚ軸は平らになっていた。
【０１６１】
１次視覚野での心磁成分除去の比較結果を図５６〜図５９に示す。推定結果で４０［ｍｍ］以上のずれは省いた。
【０１６２】
図５６は、１次視覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。Ｘ、Ｙ、Ｚ軸が１０［ｍｍ］以上ずれていた。
【０１６３】
図５７は、１次視覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。図５６と比べて、時間と共に位置が変化しなくなってきた。Ｘ、Ｙ、Ｚ軸は１０〜２０［ｍｍ］間でずれている。
【０１６４】
図５８は、１次視覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。Ｚ軸が大きく変化した。Ｘ、Ｙ軸は元のデータに近づいた。
【０１６５】
図５９は、１次視覚野における心磁成分を除去したデータを用いた等価電流双極子位置の差を説明する図である。加算平均は比較するデータと同じ回数行った。Ｘ、Ｚ軸は比較元のデータと比べて差がなくなってきた。Ｙ軸は大きくずれていた。
【０１６６】
ＭＲＩ座標系のＸＹＺ軸の標準偏差の平均値を求め、加算平均回数の関係をまとめた（図６０)。
【０１６７】
図６０は、加算平均回数とＭＲＩ座標系ＸＹＺ軸の標準偏差を示す図である。加算平均回数が増えると、標準偏差の変化が一定になっていた。１次聴覚野では、加算平均回数５０回以降は大きな変化が無くなった。１次視覚野では加算平均回数１００回以降は大きな変化がなくなった。
【０１６８】
５．５．考察
本手法は、心磁成分のみを除去することが可能である。その為、加算平均は、ある程度必要である。１次視覚野の反応は、１次聴覚野の反応に比べて十分に小さい（図５０、図５１）ので、加算平均回数が少ない１次視覚野の結果は、２０〜３０［ｍｍ］ずれることになったと考える。
【０１６９】
１次聴覚野の結果は、Ｚ軸が大きく変化した。これは、３．４章で示唆した結果と同じになった。１次視覚野でもＺ、Ｙ軸に影響をおよぼす結果があった（図５８、図５９)。まとめると、心磁成分は図６１のように、推定位置に影響を与えることが分かった。
【０１７０】
図６１は、心磁成分の等磁場線図が等価電流双極子推定に与える影響を説明する図である。赤は磁場の吹出、緑は吸い込みをあらわす。青枠を使用して等価電流双極子推定を行うと、特にＺ軸方向に影響を特に受ける。赤枠を使用したチャネルは、吹出と吸い込みの境目の為に、心磁成分の影響を大きく受ける。
【０１７１】
視覚、聴覚の比較結果で加算平均回数が増すと、時間と共に差が変化しなくなり、一定の差だけが現れるようになった（図５５、図５９）。また、図６０から加算平均回数を増しても、一定の差だけが生じていることが分かる。この一定の差は、心磁成分だと考えられる。
【０１７２】
解析データは、２−５０［Ｈｚ］のバンドパスフィルタをかけているので、直流のような成分は減衰されている。この事から、心磁成分は基線をゆがませる一要因になっているということがいえる。通常用いられる基線補正は、計測データを直接変化させる為に、どこを基準に減算するかによって、大きく計測データをゆがませてしまう。また、基線が揺らいでいるデータに対しては有効ではない。本手法を用いることにより、基線をゆがませている一要因の心磁の除去と全区間に対しての基線補正をかけている事と同じになる。
【０１７３】
本手法を用いると、心磁が作る基線成分と推定位置への影響を除外できる。本実験では、１００回の脳磁場加算平均データを用いることにより、表２で表す程度の影響を除外できた。
（表２）
心磁成分を除去したことによる等価電流双極子の平均位置
実験の種類 MRI X軸[mm] MRI Y軸[mm] MRI Z軸[mm]
聴覚 0.1±1.6[mm] 2.4±2.4[mm] 1.5±3.2[mm]
視覚 4.6±2.8[mm] 2.0±2.0[mm] 16.4±5.6[mm]
【０１７４】
６．まとめ
本研究により、可能となった事を下記に示す。
シールドルーム内で実験のタスクを行っても、心拍を計測できる心電計の開発、
心電を用いて加算平均することにより、心磁成分を抽出するプログラムの開発、
計測されたSQUIDチャネルデータから心磁成分を除去するプログラムの開発。
【０１７５】
また、明らかになった事を下記に示す。
MEGでも心磁成分が計測され、フィルタ処理等では除去できない。
心磁成分が与える影響は、測定部位により異なることを示した。
頭部を上から見て、左右の位置はＭＲＩ座標系Ｚ軸に影響を及ぼし、真ん中の位置は、ＭＲＩ座標系Ｚ、Ｙ軸に影響を及ぼす。
心磁成分は加算平均を行うと、基線のゆらぎの一要因となってあらわれる。
本手法を用いることで、心磁成分を除去し、その効果示すことができた。
【０１７６】
７．参考文献
[1]樋口正法、磁束計支援システム , Med.Imag.Tech.Vol.1l No.2 June 1993.
[2]入戸野宏、堀忠雄、心理学研究における事象関連電位 (ERP) の利用、広島大学総合科学部紀要IV理系編、第26巻、pp15-31.
[3]小杉幸夫、武者利光（2000）、電子情報通信工学シリーズ 生体情報工学、森北出版株式会社、第1版、pp13-14.
[4]野村亮輔(2006)、Shepard Tone聴音時の脳磁場反応、平成18 年度金沢工業大学人間情報工学科卒業論文,pp7-25.
[5]小林哲生、栗城真也、多チャネル誘発脳磁界計測における自発脳磁界成分低減の検討、信学技報1996-06、pp7-12.
[6]池田思朗、村田昇、Independent Component Analysis を用いた MEG データの解析、電子情報通信学会技術研究報告、Vol.98、No.129、pp. 29-36.
[7]甲斐島 武、水野（松本） 由子、伊達 進、篠崎 和弘、下條 真司、MEGデータの独立成分分析（ICA）を用いた信号成分のクラスタリング解析、日本生磁気学会誌 特別号、Vol. 16、 No. 1、pp. 256-257.
[8]平野三千代、運動負荷後QRS時間およびQT間隔の回復時間と心電図異常との関係、岩手女子看護短期大学紀要、Vol.2(19950831)、pp15-21.
[9]山岸茂孝、津田良一、精神作業負荷心電図の周波数解析、第19回熊本県産学官技術交流会講演論文集、2005年1月、pp206-207
[10]トランジスタ技術編集部、上田智章(2006)、トランジスタ技術 2006年1月号、株式会社 CQ出版、43(1)、2006年1月、pp223-233.
[11]山口巖 監修、塚田啓二 編著、心磁図の読み方、株式会社コロナ社、2006年、第１版第１刷、pp10-12.
[12]山口武志、ベッドサイドのBasic Cardiology 心筋細胞の電気生理学―イオンチャネルから、心電図、不整脈へ―、メディカル・サイエンス・インターナショナル、２００６年、第１版第５刷、pp23-90.
[13]吉村正治(1970)、脈波判読の実際、株式会社中外医学社、第4版、pp1-31.
[14]DALE DUBIN,M.D.、図解心電図テキスト、株式会社 文光堂、2007年、第6版第３刷、p66.
[15]DALE DUBIN,M.D.、図解心電図テキスト、株式会社 文光堂、2007年、第6版第3刷、p36-54.
[16]イブ トーマス、中村 尚五、プラクティス デジタル信号処理、東京電機大学出版局、2005年、第１版４刷、pp132-154.
[17]村 上 正 浩、鈴 木 博 之、内 藤 茂 昭(2006)、日立心臓磁気計測システム MC-6400、Lab.Clin.Pract.、24(2)、pp132-139.
[18]トランジスタ技術編集部、山口晶大(2002)、トランジスタ技術 2002年10月号、株式会社 CQ出版,39(10)、2002年10月、pp173-186.
[19]冨田教幸(2005)、音階知覚に関連する脳磁場応答、平成17 年度金沢工業大学人間情報工学科卒業論文、pp41-56.
[20]山口和夫(2003)、Cognitive Timerに関連する脳磁場応答の計測、平成15年度金沢工業大学人間情報工学科卒業論文、pp40-50.
[21]山内光哉、心理・教育のための統計法<第２版>、株式会社サイエンス社、第２版 第８刷、pp68-78.
【０１７７】
８．付録
８．１．製作したプログラムについて
プログラムの使用方法について述べる。ＳＱＵＩＤチャネルの数、全チャネル数のデフォルト値は人間情報システム研究所のＭＥＧ計測装置にあわせている。
【０１７８】
図６２は、心磁成分除去プログラムの画面説明図である。
黄線の枠のOpen Binボタン、Searchボタン、SUB Preボタン、SUB Outボタンを順に押すことにより心磁成分が除去される。赤枠は、閾値の設定や読み込む心電チャネルを設定する。
【０１７９】
使用方法
１．ECG CHとフィルタを選ぶ。
２．ファイル読み込み(Open Binボタン)を行う。
３．Searchボタンを押すと心電波形ウインドウとトリガリストが更新される。
４．検出されたトリガを選択するか、波形の位置を選択する。すべて、正常にＥＣＧ波形のＲ波成分に赤線がかかっていた場合は５に進む。
トリガリストの数が極端に少ない場合、又は極端に多い場合は下記のことを試す。
Ｉ．心電波形が反転している場合は、反転をチェックしSearchボタンを押す。
ＩＩ．心電波形の閾値（赤線）が高すぎる、低すぎる場合はECG Coeffを設定する。本プログラムの閾値は、心電波形のPeak Threshold Ratioを求めている。ECG Coefの値を大きくすると閾値は高くなる。
ＩＩＩ．心電波形がＲ波ではない所に、トリガが設定されていた場合は、そのトリガを選択しDeleteボタンを押す。すると、そのトリガは削除される。
ＩＶ．読み込んだ波形が、心電ではなかった場合は、１からやり直す。
５．Sub Preボタンを押すと、心磁成分除去の為の前処理とMEG CHで選択したデータに対して、心磁成分を除去した結果を示せるようになる。
６．Sub Outボタンで心磁成分が除去されたファイルが出力される。
読み込んだファイルの名前が”abc.bin”だった場合は、”abc_MCG.con”というファイル名で、同じフォルダに出力されている。
【０１８０】
トリガリストのインポート、エクスポート
トリガリストのImport、ExportボタンはMEG Laboratoryのトリガのインポート、エクスポートに対応している。MEG Laboratoryを用いて、トリガをインポート、エクスポートすることも可能である。
【０１８１】
本プログラムで、読み込むことが出来る脳磁場データを図６３に示す。MEG Laboratoryは、図６３のような構造を持つバイナリファイルをエクスポート可能である。
【０１８２】
図６３は、プログラムが読み込むことの出来るバイナリファイル構造を説明する図である。人間情報システム研究所のＭＥＧ計測装置は、１サンプルあたり１９２ｃｈのデータを記録している。本プログラムもこの仕様にあわせた。１ｃｈ、１サンプルのデータはdouble型の８バイトである。
【０１８３】
本プログラムでは、計測した心電波形に対して、デジタルフィルタをかけて、基線の揺らぎとハムノイズを除去する。ＩＩＲフィルタでバンドパスフィルタを行うと、位相ひずみの為に、波形がひずんでしまう（図６４）。その為に、本プログラムでは、ＦＩＲ型のバンドパスフィルタとＩＩＲ型のバンドエリミネーションフィルタの２種類のフィルタを実装している。ＩＩＲフィルタのフィルタ係数の算出には、株式会社デジタルフィルターのDSPLinksを用いた。ＦＩＲフィルタのタップ係数の算出は、参考文献[14]を参考にして行った。
【０１８４】
図６４は、位相ひずみによる波形ゆがみを説明する図である。Raw Dataに対して、３−５０［Ｈｚ］のバンドパスフィルタをかける。ＩＩＲフィルタは、ゆがみ（図の※）が生じてしまっている。ＦＩＲフィルタは、ゆがみが生じない。
【０１８５】
本プログラムは、心電と同時記録することにより、心磁成分を除去することが出来る。その為に、心電を同時計測していないデータには適用できない。そこで、変化がはっきりとわかる、チャネルから、心電のトリガ抽出を試みた。波形の変化が顕著なチャネル（１０８ＣＨ）と、典型的な心磁波形からピアスンの積率相関係数[21]を求め、トリガ抽出を試みた。
【０１８６】
典型的な心磁波形は、被験者４人分の平均値を用いた（図６５）。変化の顕著なＱＲＳに相当する部分を用いた。
【０１８７】
図６５は、典型的な心磁波形を示す図である。計測した被験者４人分の心磁波形の平均値を求めて、変化の顕著な部分のみを用いた。
【０１８８】
図６６は、相関を用いてトリガの抽出を試みた結果を説明する図である。青線が心電トリガの時刻を表している。赤線は、１０８ｃｈと図６５の波形との相関を表す。ｒ＝０．８５以上の時をプロットした。
【０１８９】
図６６の結果から、心電より求めたトリガの時間と、相関より求めた時間が一致している部分もあるが、一致しない部分も出力された。単純に相関を求めるだけでは、トリガ抽出は出来なかった。何らかの拘束条件を用いて、探索する必要があると考える。例えば、ヒトの心電間隔（Ｒ−Ｒ間隔）は、安静時ほぼ一定になる特性を利用し除外するとよいかもしれない。他にも変化のはっきりとわかる他のチャネルを複数用いれば、ＳＱＵＩＤセンサから計測されるデータだけで、心電トリガの抽出が可能になるかもしれない。
【０１９０】
８．２．本手法を用いたハムノイズの除去
ハムノイズの原因は交流電源である。そこで、交流電源の周期を手に入れれば、本手法によりSQUIDセンサに混入するハムノイズのみ除去が出来るはずである。計測時に６０［Ｈｚ］のバンドエリミネーションフィルタをかけていない磁場データの周波数解析結果を示す（図６７）。
【０１９１】
図６７は、計測時に６０［Ｈｚ］のバンドエリミネーションフィルタをかけていない磁場データの周波数分布を示す図である。解析に使用したＳＱＵＩＤセンサ４ｃｈ分の平均した周波数分布を示す。６０［Ｈｚ］の成分とその高調波である１２０［Ｈｚ］、１８０［Ｈｚ］にピークが見られる(図の※)。
【０１９２】
オフラインにより、６０［Ｈｚ］のバンドエリミネーションフィルタをかけてハムノイズを除去した結果を図６８に、本手法を用いてハムノイズを除去した結果を図６９に示す。
【０１９３】
図６８は、バンドエリミネーションフィルタを用いたハムノイズ除去の結果を説明する図である。６０［Ｈｚ］成分のみが減衰している（赤矢印）。しかし、高調波成分は除去されない。
【０１９４】
図６９は、本手法を用いたハムノイズ除去の結果を説明する図である。６０［Ｈｚ］、高調波成分の１２０［Ｈｚ］、１８０［Ｈｚ］のピークが消えた。バンドエリミネーションフィルタを用いた結果と比較して、６０［Ｈｚ］の成分のみ減衰していない。
【０１９５】
ハムノイズの原因である交流電源と位相を完全に固定して除去している為に、６０［Ｈｚ］の脳磁場成分は減衰させることなく除去できる。その為に、ハムノイズ除去にも本手法が有効だと考えられる。
【０１９６】
なお、プログラムの開発にはMicrosoft CorporationのMicrosoft Visual C++ 2005を用いた。聴覚実験用の音刺激作成プログラムは、参考文献[19]のプログラムを使用した。
【０１９７】
以上、本発明を実施の形態をもとに説明した。実施の形態は例示であり、それらの各構成要素や各処理プロセスの組合せにいろいろな変形例が可能なこと、またそうした変形例も本発明の範囲にあることは当業者に理解されるところである。
【図面の簡単な説明】
【０１９８】
【図１】本実施の形態に係る脳磁場計測装置の構成図である。
【図２】ＭＥＧで計測される信号を説明する図である。
【図３】ＭＥＧで計測された生体現象によるノイズを説明する図である。
【図４】ＭＥＧで計測された環境が起因となっているノイズを説明する図である。
【図５】ランダムなノイズの自己相関を説明する図である。
【図６】時間相関のあるノイズを説明する図である。
【図７】MEG実験中の風景を説明する図である。
【図８】本研究が提案する心磁成分を除去する手法を説明する図である。
【図９】基線補正を説明する図である。
【図１０】心臓の各部位の名前と特徴を説明する図である。
【図１１】実験タスクの例を説明する図である。
【図１２】脈波観測装置を説明する図である。
【図１３】脈波観測装置の回路図である。
【図１４】製作した脈波観測装置を説明する図である。
【図１５】標準１２誘導の電極貼り付け位置を説明する図である。
【図１６】典型的な心電波形を説明する図である。
【図１７】製作した心電計の回路図である。
【図１８】製作した心電計を説明する図である。
【図１９】製作した電極リードと日本光電製の電極リードを説明する図である。
【図２０】空気動圧センサを用いた脈波計測の結果を説明する図である。
【図２１】光センサを用いた脈波計測の結果を説明する図である。
【図２２】フィルタをかける前の心電波形とバンドパスフィルタ後の心電波形を説明する図である。
【図２３】胸部誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。
【図２４】肢誘導で電極貼り付け位置の違いによる波形の違いを説明する図である。
【図２５】心電計測中に電気刺激、ボタン押し動作をしたときの心電波形を説明する図である。
【図２６】光センサを用いた脈波観測装置と製作した心電計の心拍計測を説明する図である。
【図２７】電極の貼り付け位置を説明する図である。
【図２８】全頭型１６０チャンネル脳磁計を説明する図である。
【図２９】加算平均によって抽出した心磁波形と等磁場線図を説明する図である。
【図３０】各ＳＱＵＩＤセンサに含まれる心磁成分を説明する図である。
【図３１】心電波形と計測磁場データとの比較を示す図である。
【図３２】心電波形と計測した心磁図の比較を示す図である。
【図３３】被験者３人の心磁成分抽出結果を説明する図である。
【図３４】心磁波形の相似を説明する図である。
【図３５】心磁波形のパワースペクトル分布を説明する図である。
【図３６】抽出成分に対して等価電流双極子推定をおこなった結果を説明する図である。
【図３７】心磁成分と特定部位の等価電流双極子との混合波形を説明する図である。
【図３８】１次聴覚野、１次視覚野に電流双極子があった時の等磁場線図を説明する図である。
【図３９】解析に用いたチャネルを説明する図である。
【図４０】心磁成分が等価電流双極子に与える影響を説明する図である。
【図４１】心磁成分が等価電流双極子に与える影響を説明する図である。
【図４２】心磁成分が１次聴覚野近傍、１次視覚野近傍に与える等磁場線図である。
【図４３】心磁成分を除去するプログラムのスクリーンショットを示す図である。
【図４４】プログラムのフローチャートである。
【図４５】本プログラムを用いた心磁成分の除去結果（プログラムのスクリーンショット）を説明する図である。
【図４６】呈示した音刺激を説明する図である。
【図４７】呈示した視覚刺激を説明する図である。
【図４８】ＭＲＩ座標系を説明する図である。
【図４９】心磁成分除去の評価方法を説明する図である。
【図５０】聴覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。
【図５１】視覚刺激によって誘発された磁場応答と等価電流双極子位置を説明する図である。
【図５２】１次聴覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。
【図５３】１次聴覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。
【図５４】１次聴覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。
【図５５】１次聴覚野における等価電流双極子位置の差を説明する図である。
【図５６】１次視覚野における加算平均回数２０回の場合の等価電流双極子位置の差を説明する図である。
【図５７】１次視覚野における加算平均回数５０回の場合の等価電流双極子位置の差を説明する図である。
【図５８】１次視覚野における加算平均回数１００回の場合の等価電流双極子位置の差を説明する図である。
【図５９】１次視覚野における心磁成分を除去したデータを用いた等価電流双極子位置の差を説明する図である。
【図６０】加算平均回数とＭＲＩ座標系ＸＹＺ軸の標準偏差を示す図である。
【図６１】心磁成分の等磁場線図が等価電流双極子推定に与える影響を説明する図である。
【図６２】心磁成分除去プログラムの画面説明図である。
【図６３】プログラムが読み込むことの出来るバイナリファイル構造を説明する図である。
【図６４】位相ひずみによる波形ゆがみを説明する図である。
【図６５】典型的な心磁波形を示す図である。
【図６６】相関を用いてトリガの抽出を試みた結果を説明する図である。
【図６７】計測時に６０Ｈｚのバンドエリミネーションフィルタをかけていない磁場データの周波数分布を示す図である。
【図６８】バンドエリミネーションフィルタを用いたハムノイズ除去の結果を説明する図である。
【図６９】本手法を用いたハムノイズ除去の結果を説明する図である。
【符号の説明】
【０１９９】
１０ 生体信号計測部、 ２０ ノイズ信号抽出部、 ３０ ノイズ発生時刻取得部、 ４０ ノイズ除去処理部、 ５０ 生体計測信号記憶部、 ６０ ノイズ信号記憶部、 １００ 脳磁場計測装置。
【特許請求の範囲】
【請求項１】
生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、
前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、
メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとを含むことを特徴とするノイズ除去方法。
【請求項２】
前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取る際、当該ノイズ信号の時間方向の揺らぎに応じて加算平均を取る時間の範囲を一部重複させて加算平均を取ることを特徴とする請求項１に記載のノイズ除去方法。
【請求項３】
生体から発生する信号を計測する計測部と、
前記計測部による生体計測信号の計測の際に発生する時間相関性のあるノイズ信号の発生時刻に関する情報を取得する取得部と、
前記計測部により計測された生体計測信号を前記時間相関性のあるノイズ信号の発生時刻に関する情報とともに格納する記憶部と、
前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するノイズ信号抽出部と、
前記時間相関性のあるノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得するノイズ除去処理部とを含むことを特徴とする生体情報計測装置。
【請求項４】
生体から発生する脳磁場信号を計測する計測部と、
前記脳磁場信号を計測する際の心拍の発生時刻に関する情報を取得する取得部と、
前記計測部により計測された脳磁場信号を前記心拍の発生時刻に関する情報とともに格納する記憶部と、
前記脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出するノイズ信号抽出部と、
心拍の発生時刻に合わせて、抽出された前記ノイズ信号を前記脳磁場信号から減算することにより、前記ノイズ信号が除去された脳磁場信号を取得するノイズ除去処理部とを含むことを特徴とする脳磁場計測装置。
【請求項５】
生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、
前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、
メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとをコンピュータに実行させることを特徴とするプログラム。
【請求項１】
生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、
前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、
メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとを含むことを特徴とするノイズ除去方法。
【請求項２】
前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取る際、当該ノイズ信号の時間方向の揺らぎに応じて加算平均を取る時間の範囲を一部重複させて加算平均を取ることを特徴とする請求項１に記載のノイズ除去方法。
【請求項３】
生体から発生する信号を計測する計測部と、
前記計測部による生体計測信号の計測の際に発生する時間相関性のあるノイズ信号の発生時刻に関する情報を取得する取得部と、
前記計測部により計測された生体計測信号を前記時間相関性のあるノイズ信号の発生時刻に関する情報とともに格納する記憶部と、
前記生体計測信号に対して、前記時間相関性のあるノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するノイズ信号抽出部と、
前記時間相関性のあるノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得するノイズ除去処理部とを含むことを特徴とする生体情報計測装置。
【請求項４】
生体から発生する脳磁場信号を計測する計測部と、
前記脳磁場信号を計測する際の心拍の発生時刻に関する情報を取得する取得部と、
前記計測部により計測された脳磁場信号を前記心拍の発生時刻に関する情報とともに格納する記憶部と、
前記脳磁場信号に対して、心拍の発生時刻に同期させて加算平均を取ることにより、心磁によるノイズ信号を抽出するノイズ信号抽出部と、
心拍の発生時刻に合わせて、抽出された前記ノイズ信号を前記脳磁場信号から減算することにより、前記ノイズ信号が除去された脳磁場信号を取得するノイズ除去処理部とを含むことを特徴とする脳磁場計測装置。
【請求項５】
生体から計測される生体計測信号を時間相関性のあるノイズ信号の発生時刻に関する情報とともにメモリに記録するステップと、
前記生体計測信号に対して、メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に同期させて加算平均を取ることにより、前記時間相関性のあるノイズ信号を抽出するステップと、
メモリに記録された前記時間相関性のあるノイズ信号の発生時刻に関する情報を参照して当該ノイズ信号の発生時刻に合わせて、抽出された前記ノイズ信号を前記生体計測信号から減算することにより、前記ノイズ信号が除去された生体計測信号を取得し、メモリに記録するステップとをコンピュータに実行させることを特徴とするプログラム。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【図１９】
【図２０】
【図２１】
【図２２】
【図２３】
【図２４】
【図２５】
【図２６】
【図２７】
【図２８】
【図２９】
【図３０】
【図３１】
【図３２】
【図３３】
【図３４】
【図３５】
【図３６】
【図３７】
【図３８】
【図３９】
【図４０】
【図４１】
【図４２】
【図４３】
【図４４】
【図４５】
【図４６】
【図４７】
【図４８】
【図４９】
【図５０】
【図５１】
【図５２】
【図５３】
【図５４】
【図５５】
【図５６】
【図５７】
【図５８】
【図５９】
【図６０】
【図６１】
【図６２】
【図６３】
【図６４】
【図６５】
【図６６】
【図６７】
【図６８】
【図６９】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【図１９】
【図２０】
【図２１】
【図２２】
【図２３】
【図２４】
【図２５】
【図２６】
【図２７】
【図２８】
【図２９】
【図３０】
【図３１】
【図３２】
【図３３】
【図３４】
【図３５】
【図３６】
【図３７】
【図３８】
【図３９】
【図４０】
【図４１】
【図４２】
【図４３】
【図４４】
【図４５】
【図４６】
【図４７】
【図４８】
【図４９】
【図５０】
【図５１】
【図５２】
【図５３】
【図５４】
【図５５】
【図５６】
【図５７】
【図５８】
【図５９】
【図６０】
【図６１】
【図６２】
【図６３】
【図６４】
【図６５】
【図６６】
【図６７】
【図６８】
【図６９】
【公開番号】特開２００９−１９５５７１（Ｐ２００９−１９５５７１Ａ）
【公開日】平成２１年９月３日（２００９．９．３）
【国際特許分類】
【出願番号】特願２００８−４２１１９（Ｐ２００８−４２１１９）
【出願日】平成２０年２月２２日（２００８．２．２２）
【出願人】（５９３１６５４８７）学校法人金沢工業大学 (202)
【Ｆターム（参考）】
【公開日】平成２１年９月３日（２００９．９．３）
【国際特許分類】
【出願日】平成２０年２月２２日（２００８．２．２２）
【出願人】（５９３１６５４８７）学校法人金沢工業大学 (202)
【Ｆターム（参考）】
[ Back to top ]