JP2004041670A - 生体機能情報の分別抽出方法およびその装置 - Google Patents
生体機能情報の分別抽出方法およびその装置 Download PDFInfo
- Publication number
- JP2004041670A JP2004041670A JP2002265051A JP2002265051A JP2004041670A JP 2004041670 A JP2004041670 A JP 2004041670A JP 2002265051 A JP2002265051 A JP 2002265051A JP 2002265051 A JP2002265051 A JP 2002265051A JP 2004041670 A JP2004041670 A JP 2004041670A
- Authority
- JP
- Japan
- Prior art keywords
- biological function
- signal
- function information
- measurement data
- physiological
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Images
Landscapes
- Magnetic Resonance Imaging Apparatus (AREA)
- Measuring And Recording Apparatus For Diagnosis (AREA)
Abstract
【課題】局所神経活動に伴う信号変化とその他の起源を有する生理学的変動が混在する場合においても、信号を生理学的起源の異なる構成成分に精度よく分別し、構成成分別の情報を提供する。
【解決手段】生体機能データを計測し(1)、必要に応じて解析領域の選択(2)、前処理(3)および入力信号の生成(4)を行い、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築する(5)。これにより生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する。抽出結果は、例えば表示装置に表示される。
【選択図】 図1
【解決手段】生体機能データを計測し(1)、必要に応じて解析領域の選択(2)、前処理(3)および入力信号の生成(4)を行い、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築する(5)。これにより生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する。抽出結果は、例えば表示装置に表示される。
【選択図】 図1
Description
【0001】
【発明の属する技術分野】
本発明は、磁気共鳴映像法や近赤外線分光法などの生体機能計測データから、生理学的起源の異なる生体機能に由来する構成成分を分別抽出し、各成分の情報を表示する生体機能情報の分別抽出方法およびその装置に関する。
【0002】
【従来の技術】
磁気共鳴映像法、近赤外線分光法などの生体機能測定法は非侵襲的な脳機能計測法として近年、広く利用されている。これらの計測データには局所神経活動に伴う信号変化のほか、全身循環、呼吸、自律神経活動、血管運動などさまざまな生理学的現象に由来する信号変化の影響が混在している。生理学的現象のうち計測が容易な心拍、呼吸の影響については、これらの周波数成分と生体機能計測データに含まれる周波数成分の比較が行われ、生体機能計測データにこれらの生理学的現象の成分が含まれていることが明らかにされている(非特許文献1、2、3)。
【非特許文献1】
Ogawa S, et al Visualization of Information Processing in the Human Brain: Recent Advances in MEG and Functional MRI (EEG Suppl. 47) 1996, pp5−14.
【非特許文献2】
Mitra PP, et al, Magn Reson Med 1997;37:511−8.
【非特許文献3】
Lowe MJ, et al, Neuroimage 1998;7:119−32.
脳機能検査法としての応用では局所神経活動とこれら生理学的現象に伴う信号変化の周波数成分の違いを利用して、トレンド除去、ローパスフィルターの適用などの前処理が行われてきた(非特許文献4、5)。
【非特許文献4】
Friston KJ, et al, Human Brain Mapping 1995;2;189−210.
【非特許文献5】
Strupp JP, Neuroimage 1996;3;S607.
従来は局所神経活動に伴う変化を抽出するためには神経活動賦活時の計測値と基準となる計測値の間の統計的検定を行い有意差のある部分を抽出する方法や、典型的な信号変化パターンとの相関を求める方法などが用いられてきた(非特許文献5)。さらにこれらの手法を拡張した一般線形モデルも用いられている(非特許文献4)。また、関心領域の反応パターンを推定するためには賦活パラダイムに同期した加算平均や事前に予測された変化パターンに基づくフィッティングなどが用いられてきた。最近では独立成分分析を用いて由来の異なる信号成分の分別が試みられている(非特許文献6)。
【非特許文献6】
Gu H, et al, Neuroimage 2001;14:1432−43.
一方、その他の生理学的現象に由来する信号は外乱として認識されており、これらの生体機能計測データからの局所生理機能抽出の応用は行われていなかった。
【0003】
【発明が解決しようとする課題】
生体機能計測データに影響を与える生理学的現象としては循環、呼吸、自律神経活動、血管運動などがあげられる。前記で述べた方法では、局所神経活動に伴う信号変化とは大きく周波数が異なる生理学的変動、たとえば心拍、呼吸の影響や数分単位の血圧、動脈血二酸化炭素分圧の変動などの影響は前処理により容易に除去可能であるが、周波数帯域の重なる変化、たとえば圧受容体反射に伴う血圧変動の影響や血管運動の影響は除去できない。したがって、これらの影響を除去するにはデータサンプル数を増やして統計パワーを増加させる必要がある。しかしながらこれらの生理学的変動は賦活パラダイムとは独立した現象とは限らず、データサンプル数の増加のみではバイアスのない検定、推定結果がえられるとは限らない。また、これらの影響が大きい場合には誤った結果を導き出す可能性がある。独立成分分析を用いた場合、分離された信号成分の生理学的な解釈は必ずしも容易ではない。そして、これら生理学的現象が生体機能計測データに及ぼす影響の強度は空間的に一様ではない。また、時間的にも同期的、非同期的な影響が混在している。
【0004】
本発明は、局所神経活動に伴う信号変化とその他の起源を有する生理学的変動が混在する場合においても、信号を生理学的起源の異なる構成成分に精度よく分別し、構成成分別の情報を提供することを目的とする。また本発明は、生体機能計測データに含まれるこれら生理学的起源を異にする信号成分の時間的あるいは空間的特徴を容易に理解可能な形として提供することを目的とする。
【0005】
【課題を解決するための手段】
上記目的は、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する生体機能情報の分別抽出方法により、達成される。
ここで、抽出される生体機能は中枢神経系の機能または局所生理機能とすることができる。また、数式モデルの構築において、生体機能計測データの代表信号値を使用すること、全身生理機能の測定値を使用すること、または確率的雑音を使用することができる。
また、分別抽出された生体機能情報を表示することができる。生体機能情報の表示において、モデルに基づくシミュレーションを使用すること、ノイズ寄与率(パワー寄与率)の分布または空間的情報の表示を使用すること、モデル固有値の分布または空間的情報の表示を使用すること、または確率的雑音の強度の分布または空間的情報の表示を使用することができる。
本発明に係る生体機能情報の分別抽出装置は、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する手段を備えて構成される。本分別抽出装置は、さらに、分別抽出された生体機能情報を表示する手段を備えることができる。
【0006】
本発明では、解析に用いる生体機能計測データとしては、磁気共鳴信号、近赤外線スペクトルなどを使用する。これら生体機能計測データにはさまざまな起源を有する生理学的変動が混在している。一般に局所神経活動に伴う変化は賦活パラダイムに対応した応答を示すが、その応答は空間的、時間的に一様なものではない。その他の起源を有する生理学的変動の由来としては心拍、呼吸といった生体の有する基本的な周期的活動、自律神経活動に伴う血圧変動、動脈血二酸化炭素分圧の変動などがあり、これらについては前記生体機能計測データに対して同期的な信号変化をもたらすが、これらによる影響の強度は空間的に一様なものではない。これに対して血管運動、すなわち小血管の周期的収縮・拡張にともなう信号変化は通常、局所ごとに異なる位相を示し、非同期的である。またその強度も空間的に一様ではない。各画素またはチャンネルに対応して入力−出力関係を記述する数式モデルを構築すれば、各種生理学的変動に対する空間的、時間的応答を推定することができる。本発明では、生体機能計測によりえられた画素またはチャンネル毎の計測値を出力とし、機能賦活パラダイム、同時に計測された生理機能の計測値、生体機能計測データから作成された参照信号、確率的雑音などを入力とする数式モデルを構築することで、各入力が生体機能計測値に及ぼす影響を分別抽出する。えられた数式モデルのモデルパラメーターの表示に加え、各入力を独立に変化させたときの出力のシミュレーションを表示することで入力―出力間の時間的、空間的関係を視覚的に容易に理解可能な形で提供し、前記目的を達成する。入力成分間に相互作用がある系においても各成分が他の成分の信号変化に寄与する割合を表示することで各成分間のフィードバック関係を理解容易な形で提供する。また、明らかな入力を特定しえない生理学的変動、たとえば血管運動に起因する信号変化に対しては確率的雑音を入力としてモデルを構築し、そのモデルパラメーターの表示に加え、その生理学的変動成分の特徴を理解容易な形にして提供する。また、本発明では連続した測定データをモデル推定用とモデル評価用に分割して使用することで、構築されたモデルの妥当性を検証することができる。あるいは別個に収集されたデータを結合してモデルを構築することもできるので、データ収集に関する制約が少なくなる。
【0007】
【発明の実施の形態】
以下、本発明の実施の形態を例とともに説明する。図1は本発明の実施手順を示すフローチャートである。本フローチャートの各ステップ1〜6は次のような処理が行われる。
【0008】
(1)生体機能データの計測
図2(a)、(b)に生体機能および生理学的データの収集方法を示す。磁気共鳴画像(図中a1)では高い時間分解能がえられるエコープラナー法を使用することが望ましい。近赤外線スペクトル(図中a2)では空間的情報をえるため、多チャンネル計測を使用することが望ましい。脳機能計測では被験者に対し刺激ないし課題遂行を行う機能賦活期間と刺激ないし課題を与えない対照期間のデータを収集する。図中、右方向は時間軸であり、賦活期間と対照期間が交互に付与される。本発明ではデータ収集に際して従来の方法のように一定の機能賦活期間と対照期間の組み合わせを反復して行う必要はなく、上記期間の長さは一回のデータ収集期間の間であっても、図2に示すように、試行毎に変えて差し支えない。局所生理機能の情報抽出のみを目的とする場合には、刺激ないし課題を与える必要はない。生体機能計測データの収集と平行して、生理学的変動の参照データとして脈波、連続血圧(図中b1)、胸腹壁運動、呼気中二酸化炭素濃度(図中b2)などを計測し、生理学的データの収集を行う。なお、あとで示すように生体機能計測データ計測値の代表値を参照データとして使用する場合は、これらの生理学的データの収集は必要とはしない。
【0009】
(2)解析領域の選択
磁気共鳴画像を生体機能計測データとして使用する場合は解析を行う領域を脳が存在する部分に制限する。その際には、脳が存在する部分と背景の領域の磁気共鳴信号の強度差を利用し、当該信号強度に閾値を設定することで両者を判別することができる。なお、データによってはこの処理を省略しても差し支えない。
近赤外線スペクトルを生体機能計測データとして使用する場合にはこの処理は必要としない。
【0010】
(3)前処理
(1)、(2)の方法で生成したデータに対し、直流および超低周波数成分を除去する平均値除去、トレンド除去処理を画素またはチャンネル毎に行う。また、高周波数成分を除去するローパスフィルター処理を行う。さらに関心のある周波数帯域に合わせてダウンサンプリングを行う。なお、データによってはこれらの処理を省略しても差し支えない。
【0011】
(4)入力信号の生成
入力信号として使用する脳機能賦活パラダイムデータと生理学的変動の参照データを生成する。脳機能賦活パラダイムデータは機能賦活期間と対照期間それぞれに任意の値を割り付けて生成するが、期間全体の平均値が0となるように値を設定することが望ましい。生理学的変動の参照データとしては連続血圧、呼気中二酸化炭素濃度などの実測値を使用するが、データ収集時の被験者の負担が少なく、後に述べる数式モデルの構造を簡略化するための方法として(1)、(2)、(3)の処理をへた生体機能計測データの代表値を使用することも可能である。本実施例では生体機能計測データの代表値として全脳平均値を使用する。なお、病態により正常と異なる信号変化を示す部分が生体機能計測データに含まれている場合には、当該の画素またはチャンネルの計測値を除外して処理を行う。多断面の磁気共鳴画像を用いる場合には断面ごとの撮像タイミングのずれを補正するため、断面毎の撮像タイミングのずれを調整した補間処理を行った後、平均値を算出する。これら測定値に対し、(3)と同様にして直流および超低周波数成分を除去する平均値除去、トレンド除去処理を行う。また、高周波数成分を除去するローパスフィルター処理を行う。なお、データによってはこの処理を省略しても差し支えない。生理学的変動の参照データに遅延がある場合、たとえばサイドストリームとして採取した呼気ガスサンプルを用いて測定した二酸化炭素濃度を使用する場合などには、あらかじめ遅延時間を補正した計測値を入力信号として使用する。生理学的変動の参照データのサンプリングレートが生体機能計測データのそれと異なる場合には両者が同じサンプリングレートとなるようにリサンプリング処理により調節する。
【0012】
(5)数式モデルの生成
利用する数式モデルとしては線形、時不変性モデルに限らず、非線形、時変性モデルを含むあらゆるブラックボックスモデルを使用することができる。生体の活動は一般的には非線形、時変性であるが、安定した生理学的状態、比較的短い時間の範囲では線形、時不変性システムと仮定して取り扱っても大きな誤差は生じない。本実施例ではモデルパラメーターの計算が容易な自己回帰システムモデル(auto−regressive exogenous model; ARX model)を使用して説明する(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著.MATLABによる制御のためのシステム同定,東京電機大学出版局.Ljung L, System Identification Toolbox User’s Guide Version 5, The Mathworks, Inc)。数式を次に示す。
【0013】
【数1】
【0014】
ただし、y(k)はシステムの出力信号、u(k)はシステムの入力信号、nkはむだ時間、w(k)は平均値0、有限な分散を持つ白色雑音、A(q)とB(q)はシフトオペレーターqの多項式である。
画素またはチャンネル毎の計測値を出力、脳機能賦活パラダイムデータと生理学的変動の参照データを入力としてモデルパラメーターの推定を行う。ARXモデルパラメーターの推定には最小二乗推定法を使用するが、これ以外のモデル、たとえば自己回帰移動平均システムモデル(auto−regressive moving average exogenous model; ARMAX model)を使用する場合にはモデルパラメーター推定法として最尤推定法などを使用する(中溝高好著,信号解析とシステム同定,コロナ社)。
【0015】
示適モデル次数はクロスバリデーション、最小実現、情報規範などを用いて選定する(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著,MATLABによる制御のためのシステム同定,東京電機大学出版局)。本実施例では多数の画素またはチャンネルの計測値を処理するにあたって、オペレーターの関与を少なくする方法として情報規範を用いた自動的選定、あるいは自動的に選定されたモデル次数の結果を参考にしたオペレーター選定によった。自動選定に用いる情報規範としては赤池情報量基準などを使用する(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。同様にしてむだ時間nkについても相関法を用いたインパルス応答あるいはクロスバリデーション、情報規範などで選定を行う。なお、データによってはむだ時間を0に固定しても差し支えない。
【0016】
全脳平均値を入力信号として使用した場合には、全脳平均値の変化に先行して信号変化を生ずる脳部位が存在する可能性があるため、全脳平均値の遅延時間の調整を行う。調整量はクロスバリデーションに基づく選定、あるいは既知の通過時間に基づきオペレーターの指示による選定を行う。なお、データによってはこの処理を省略しても差し支えない。
【0017】
局所生理機能の情報抽出のみを目的とする場合には、生理学的変動の参照データを入力とするモデルを使用する。あるいは複数の生理学的変動の参照データを用い、これらの間のフィードバック関係を含めて情報抽出を行う場合は多変量時系列モデルを使用する(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。明らかな入力を特定しえない生理学的変動、たとえば血管運動に起因する信号変化の特徴を抽出する場合は他の生理学的変動の参照データとともに確率的雑音を入力としてモデル化する。あるいは他の生理学的変動の影響が除外できる場合は時系列モデルを使用する。
【0018】
(6)結果の表示
画素、チャンネル毎に推定されたモデルパラメーターは多項式表現、伝達関数表現、零点−極点−ゲイン表現などに形式を変換して表示することができる。あるいは周波数伝達関数として表示することができる。零点、極点などのモデル固有値は複素平面上での分布表示あるいは周波数帯域別の空間的分布を示すマップ表示として示すことができる。
図3はヒトの視覚刺激時の頭部磁気共鳴画像データに本発明を適用し得られたモデルに基づき、脳機能賦活パラダイムに対応するゲインの空間的分布を示すマップである。ゲインの値は画像の輝度により表示している。
図4はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルの固有値(極点)の分布を複素平面上に表示したグラフである。
図5および図6は本発明を適用し得られたモデルの固有値(極点および零点)の周波数帯域別の空間的分布を示すマップである。これらの例においては周波数帯域 0.008〜0.05 Hz、0.05〜0.10 Hz、0.10〜0.15 Hz、0.15〜0.20 Hz、0.20〜0.25 Hzにおける画素毎の固有値の数を左上から右下に向かって配置された画像によって示す。なお最下段右はブランクである。固有値の数は画像の輝度により表示している。図6では、最上段右の周波数帯域 0.05〜0.10 Hzの画像のみが現われている。
確率的雑音の強度は分散あるいは標準偏差、変動係数などに変換して、その分布あるいは空間的分布を示すマップとして表示することができる。
図7はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルの確率的雑音の強度の空間的分布を示すマップである。確率的雑音の強度は画像の輝度により表示している。この例では表示スケールは画素毎の確率的雑音の分散に1を加えた値の常用対数により示している。
モデルに基づいてシミュレーションを行う場合にはインパルス応答、ステップ応答、その他、オペレーターの指定する任意の入力波形に対する応答を表示する。応答は入力開始時を基準とする経時的変化、あるいはオペレーターの指定する経過時間後の信号値の空間的分布を示すマップとして表示することができる。
【0019】
生体機能計測の実測値と入力信号に基づいたシミュレーション結果を比較することで、推定されたモデルの妥当性を検証することができる。
図8はヒトの視覚刺激時の頭部磁気共鳴画像データに本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。この例においては8秒間の強度1のパルス状入力に対する脳内各部位の0.5秒間隔の反応を左上から右下に向かって配置された画像により示す。2段目左端の画像が入力開始時、上から4段目左から4番目の画像が入力終了時に対応する。入力に対応する信号強度の変化量は画像の輝度により表示している。
図9は後頭葉の1画素の信号変化のシミュレーション結果を示すグラフである。横軸は時間経過を入力開始時点を1とする画像番号により表示している。縦軸は信号強度の変化量を示す。参照のため入力波形を重ねて表示している。
図10は上記モデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。この例においては8秒間の強度1のパルス状入力に対する脳内各部位の0.5秒間隔の反応を左上から右下に向かって配置された画像により示す。2段目左端の画像が入力開始時、上から4段目左から4番目の画像が入力終了時に対応する。入力に対応する信号強度の変化量は画像の輝度により表示している。
図11は頭頂葉の1画素の信号変化のシミュレーション結果を示すグラフである。横軸は時間経過を入力開始時点を1とする画像番号により表示している。縦軸は信号強度の変化量を示す。参照のため入力波形を重ねて表示している。
【0020】
本実施例は、脳機能賦活パラダイムとして1種類の入力を使用する例を示したが、複数の入力を用いる脳機能賦活パラダイムに対しても適用することができる。
本実施例は、1断面の磁気共鳴画像を用いているが、多断面のデータでも各断面に対して同様の手順で取り扱うことができる。多断面の場合には断面毎の撮像タイミングのずれを補正して処理を行う。
本実施例では、全データ収集後の一括処理を例にとって説明したが、逐次計算アルゴリズムを使用することにより(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著,MATLABによる制御のためのシステム同定,東京電機大学出版局)、オンラインあるいはリアルタイム処理が可能である。
本実施例ではヒトを対象とした脳機能賦活情報の抽出を例にとって説明したが、本発明は脳以外の臓器、組織およびヒト以外の動物を対象として適用することができる。
【0021】
時系列モデルを使用した場合にはモデルパラメーターの表示または形式を変換した表示、固有値の表示、パワースペクトル、クロススペクトル、コヒーレンスなどの周波数領域の情報に変換した表示、あるいは確率的雑音の強度の表示をすることができる。また、モデルに基づくシミュレーションを表示することができる。多変量自己回帰モデルにおいて、各変量の雑音がお互いに無相関とみなせる場合には、変量間にフィードバックが存在する系であっても、閉ループの周波数応答関数、ノイズ寄与率(またはパワー寄与率)を表示することができる(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。
図12はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(脈波、胸郭運動、呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルに基づき、横軸に周波数、縦軸に磁気共鳴信号に対する各変量のノイズ(パワー)寄与率を表示したグラフである。図中、Aは呼気終末二酸化炭素濃度、Bは胸郭運動、Cは脈波、Dは自身(磁気共鳴信号)をそれぞれ示す。
また、本発明を適用し得られたモデルに基づき、周波数帯域を選択して生体機能計測データに対する各変量のノイズ(パワー)寄与率の空間的分布をマップとして表示することができる。
図13は本発明を適用し得られたモデルに基づき、周波数帯域0.98〜1.10 Hzにおける磁気共鳴信号に対する脈波のノイズ(パワー)寄与率の空間的分布を示すマップである。脈波のノイズ(パワー)寄与率は画像の輝度により表示している。
図14は本発明を適用し得られたモデルに基づき、周波数帯域0.24〜0.36 Hzにおける磁気共鳴信号に対する胸郭運動のノイズ(パワー)寄与率の空間的分布を示すマップである。胸郭運動のノイズ(パワー)寄与率は画像の輝度により表示している。
図15は本発明を適用し得られたモデルに基づき、周波数帯域0.04〜0.16 Hzにおける磁気共鳴信号に対する呼気終末二酸化炭素濃度のノイズ(パワー)寄与率の空間的分布を示すマップである。呼気終末二酸化炭素濃度のノイズ(パワー)寄与率は画像の輝度により表示している。
本実施例では変量毎にノイズ(パワー)寄与率の空間的分布を示すマップを作成したが、各変量に色を割り付けることにより1枚のカラーマップとして表示することも可能である。
また、任意の変量に対し入力を与えた場合の応答をシミュレーションとして表示することができる(尾崎統、北川源四郎編、時系列解析の方法、朝倉書店.尾崎統、北川源四郎編、時系列解析の実際I、朝倉書店)。
【0022】
本発明に係る生体機能情報の分別抽出装置としては、例えば、コンピュータを用いて生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出するようにすることができる。分別抽出された生体機能情報はコンピュータに接続された表示装置に表示することができる。
【発明の効果】
本発明では、既知の抽出法である平均値の有意差検定や相関を用いた方法とは異なり、刺激提示あるいは課題遂行に伴う生体機能計測データの信号変化のパターンを仮定する必要がない。また、本発明では統計的検定の前提となる事象の独立性を満足するための平滑化処理を必要としないので生体機能計測データの分解能を劣化させることがない。本発明では独立成分分析とは異なり、分別された信号成分の生理学的由来について推測する必要がない。本発明では、局所神経活動以外の生理学的変動が生体機能計測データに混在する場合でも精度よく局所神経活動の情報を抽出することができる。さらに、局所神経活動以外の種々の生理機能についての情報を抽出することができる。
【図面の簡単な説明】
【図1】生理機能賦活情報の分別抽出手順を示すフローチャートである。
【図2】(a)、(b)は生体機能および生理学的データの収集方法を示す図である。
【図3】本発明を適用し得られたモデルに基づき、脳機能賦活パラダイムに対応するゲインの空間的分布を示すマップである。
【図4】本発明を適用し得られたモデルの固有値(極点)の分布を複素平面上に表示したグラフである。
【図5】本発明を適用し得られたモデルの固有値(極点)の周波数帯域別の空間的分布を示すマップである。
【図6】本発明を適用し得られたモデルの固有値(零点)の周波数帯域別の空間的分布を示すマップである。
【図7】本発明を適用し得られたモデルの確率的雑音の強度の空間的分布を示すマップである。
【図8】本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。
【図9】本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの1画素の信号変化のシミュレーション結果を示すグラフである。
【図10】本発明を適用し得られたモデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。
【図11】本発明を適用し得られたモデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの1画素の信号変化のシミュレーション結果を示すグラフである。
【図12】本発明を適用し得られたモデルに基づき、磁気共鳴信号に対する各変量のノイズ(パワー)寄与率を表示したグラフである。
【図13】本発明を適用し得られたモデルに基づき、周波数帯域0.98〜1.10 Hzにおける磁気共鳴信号に対する脈波のノイズ(パワー)寄与率の空間的分布を示すマップである。
【図14】本発明を適用し得られたモデルに基づき、周波数帯域0.24〜0.36 Hzにおける磁気共鳴信号に対する胸郭運動のノイズ(パワー)寄与率の空間的分布を示すマップである。
【図15】本発明を適用し得られたモデルに基づき、周波数帯域0.04〜0.16 Hzにおける磁気共鳴信号に対する呼気終末二酸化炭素濃度のノイズ(パワー)寄与率の空間的分布を示すマップである。
【符号の説明】
1 生体機能データの計測
2 解析領域の選択
3 前処理
4 入力信号の生成
5 数式モデルの生成
6 結果の表示
【発明の属する技術分野】
本発明は、磁気共鳴映像法や近赤外線分光法などの生体機能計測データから、生理学的起源の異なる生体機能に由来する構成成分を分別抽出し、各成分の情報を表示する生体機能情報の分別抽出方法およびその装置に関する。
【0002】
【従来の技術】
磁気共鳴映像法、近赤外線分光法などの生体機能測定法は非侵襲的な脳機能計測法として近年、広く利用されている。これらの計測データには局所神経活動に伴う信号変化のほか、全身循環、呼吸、自律神経活動、血管運動などさまざまな生理学的現象に由来する信号変化の影響が混在している。生理学的現象のうち計測が容易な心拍、呼吸の影響については、これらの周波数成分と生体機能計測データに含まれる周波数成分の比較が行われ、生体機能計測データにこれらの生理学的現象の成分が含まれていることが明らかにされている(非特許文献1、2、3)。
【非特許文献1】
Ogawa S, et al Visualization of Information Processing in the Human Brain: Recent Advances in MEG and Functional MRI (EEG Suppl. 47) 1996, pp5−14.
【非特許文献2】
Mitra PP, et al, Magn Reson Med 1997;37:511−8.
【非特許文献3】
Lowe MJ, et al, Neuroimage 1998;7:119−32.
脳機能検査法としての応用では局所神経活動とこれら生理学的現象に伴う信号変化の周波数成分の違いを利用して、トレンド除去、ローパスフィルターの適用などの前処理が行われてきた(非特許文献4、5)。
【非特許文献4】
Friston KJ, et al, Human Brain Mapping 1995;2;189−210.
【非特許文献5】
Strupp JP, Neuroimage 1996;3;S607.
従来は局所神経活動に伴う変化を抽出するためには神経活動賦活時の計測値と基準となる計測値の間の統計的検定を行い有意差のある部分を抽出する方法や、典型的な信号変化パターンとの相関を求める方法などが用いられてきた(非特許文献5)。さらにこれらの手法を拡張した一般線形モデルも用いられている(非特許文献4)。また、関心領域の反応パターンを推定するためには賦活パラダイムに同期した加算平均や事前に予測された変化パターンに基づくフィッティングなどが用いられてきた。最近では独立成分分析を用いて由来の異なる信号成分の分別が試みられている(非特許文献6)。
【非特許文献6】
Gu H, et al, Neuroimage 2001;14:1432−43.
一方、その他の生理学的現象に由来する信号は外乱として認識されており、これらの生体機能計測データからの局所生理機能抽出の応用は行われていなかった。
【0003】
【発明が解決しようとする課題】
生体機能計測データに影響を与える生理学的現象としては循環、呼吸、自律神経活動、血管運動などがあげられる。前記で述べた方法では、局所神経活動に伴う信号変化とは大きく周波数が異なる生理学的変動、たとえば心拍、呼吸の影響や数分単位の血圧、動脈血二酸化炭素分圧の変動などの影響は前処理により容易に除去可能であるが、周波数帯域の重なる変化、たとえば圧受容体反射に伴う血圧変動の影響や血管運動の影響は除去できない。したがって、これらの影響を除去するにはデータサンプル数を増やして統計パワーを増加させる必要がある。しかしながらこれらの生理学的変動は賦活パラダイムとは独立した現象とは限らず、データサンプル数の増加のみではバイアスのない検定、推定結果がえられるとは限らない。また、これらの影響が大きい場合には誤った結果を導き出す可能性がある。独立成分分析を用いた場合、分離された信号成分の生理学的な解釈は必ずしも容易ではない。そして、これら生理学的現象が生体機能計測データに及ぼす影響の強度は空間的に一様ではない。また、時間的にも同期的、非同期的な影響が混在している。
【0004】
本発明は、局所神経活動に伴う信号変化とその他の起源を有する生理学的変動が混在する場合においても、信号を生理学的起源の異なる構成成分に精度よく分別し、構成成分別の情報を提供することを目的とする。また本発明は、生体機能計測データに含まれるこれら生理学的起源を異にする信号成分の時間的あるいは空間的特徴を容易に理解可能な形として提供することを目的とする。
【0005】
【課題を解決するための手段】
上記目的は、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する生体機能情報の分別抽出方法により、達成される。
ここで、抽出される生体機能は中枢神経系の機能または局所生理機能とすることができる。また、数式モデルの構築において、生体機能計測データの代表信号値を使用すること、全身生理機能の測定値を使用すること、または確率的雑音を使用することができる。
また、分別抽出された生体機能情報を表示することができる。生体機能情報の表示において、モデルに基づくシミュレーションを使用すること、ノイズ寄与率(パワー寄与率)の分布または空間的情報の表示を使用すること、モデル固有値の分布または空間的情報の表示を使用すること、または確率的雑音の強度の分布または空間的情報の表示を使用することができる。
本発明に係る生体機能情報の分別抽出装置は、生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する手段を備えて構成される。本分別抽出装置は、さらに、分別抽出された生体機能情報を表示する手段を備えることができる。
【0006】
本発明では、解析に用いる生体機能計測データとしては、磁気共鳴信号、近赤外線スペクトルなどを使用する。これら生体機能計測データにはさまざまな起源を有する生理学的変動が混在している。一般に局所神経活動に伴う変化は賦活パラダイムに対応した応答を示すが、その応答は空間的、時間的に一様なものではない。その他の起源を有する生理学的変動の由来としては心拍、呼吸といった生体の有する基本的な周期的活動、自律神経活動に伴う血圧変動、動脈血二酸化炭素分圧の変動などがあり、これらについては前記生体機能計測データに対して同期的な信号変化をもたらすが、これらによる影響の強度は空間的に一様なものではない。これに対して血管運動、すなわち小血管の周期的収縮・拡張にともなう信号変化は通常、局所ごとに異なる位相を示し、非同期的である。またその強度も空間的に一様ではない。各画素またはチャンネルに対応して入力−出力関係を記述する数式モデルを構築すれば、各種生理学的変動に対する空間的、時間的応答を推定することができる。本発明では、生体機能計測によりえられた画素またはチャンネル毎の計測値を出力とし、機能賦活パラダイム、同時に計測された生理機能の計測値、生体機能計測データから作成された参照信号、確率的雑音などを入力とする数式モデルを構築することで、各入力が生体機能計測値に及ぼす影響を分別抽出する。えられた数式モデルのモデルパラメーターの表示に加え、各入力を独立に変化させたときの出力のシミュレーションを表示することで入力―出力間の時間的、空間的関係を視覚的に容易に理解可能な形で提供し、前記目的を達成する。入力成分間に相互作用がある系においても各成分が他の成分の信号変化に寄与する割合を表示することで各成分間のフィードバック関係を理解容易な形で提供する。また、明らかな入力を特定しえない生理学的変動、たとえば血管運動に起因する信号変化に対しては確率的雑音を入力としてモデルを構築し、そのモデルパラメーターの表示に加え、その生理学的変動成分の特徴を理解容易な形にして提供する。また、本発明では連続した測定データをモデル推定用とモデル評価用に分割して使用することで、構築されたモデルの妥当性を検証することができる。あるいは別個に収集されたデータを結合してモデルを構築することもできるので、データ収集に関する制約が少なくなる。
【0007】
【発明の実施の形態】
以下、本発明の実施の形態を例とともに説明する。図1は本発明の実施手順を示すフローチャートである。本フローチャートの各ステップ1〜6は次のような処理が行われる。
【0008】
(1)生体機能データの計測
図2(a)、(b)に生体機能および生理学的データの収集方法を示す。磁気共鳴画像(図中a1)では高い時間分解能がえられるエコープラナー法を使用することが望ましい。近赤外線スペクトル(図中a2)では空間的情報をえるため、多チャンネル計測を使用することが望ましい。脳機能計測では被験者に対し刺激ないし課題遂行を行う機能賦活期間と刺激ないし課題を与えない対照期間のデータを収集する。図中、右方向は時間軸であり、賦活期間と対照期間が交互に付与される。本発明ではデータ収集に際して従来の方法のように一定の機能賦活期間と対照期間の組み合わせを反復して行う必要はなく、上記期間の長さは一回のデータ収集期間の間であっても、図2に示すように、試行毎に変えて差し支えない。局所生理機能の情報抽出のみを目的とする場合には、刺激ないし課題を与える必要はない。生体機能計測データの収集と平行して、生理学的変動の参照データとして脈波、連続血圧(図中b1)、胸腹壁運動、呼気中二酸化炭素濃度(図中b2)などを計測し、生理学的データの収集を行う。なお、あとで示すように生体機能計測データ計測値の代表値を参照データとして使用する場合は、これらの生理学的データの収集は必要とはしない。
【0009】
(2)解析領域の選択
磁気共鳴画像を生体機能計測データとして使用する場合は解析を行う領域を脳が存在する部分に制限する。その際には、脳が存在する部分と背景の領域の磁気共鳴信号の強度差を利用し、当該信号強度に閾値を設定することで両者を判別することができる。なお、データによってはこの処理を省略しても差し支えない。
近赤外線スペクトルを生体機能計測データとして使用する場合にはこの処理は必要としない。
【0010】
(3)前処理
(1)、(2)の方法で生成したデータに対し、直流および超低周波数成分を除去する平均値除去、トレンド除去処理を画素またはチャンネル毎に行う。また、高周波数成分を除去するローパスフィルター処理を行う。さらに関心のある周波数帯域に合わせてダウンサンプリングを行う。なお、データによってはこれらの処理を省略しても差し支えない。
【0011】
(4)入力信号の生成
入力信号として使用する脳機能賦活パラダイムデータと生理学的変動の参照データを生成する。脳機能賦活パラダイムデータは機能賦活期間と対照期間それぞれに任意の値を割り付けて生成するが、期間全体の平均値が0となるように値を設定することが望ましい。生理学的変動の参照データとしては連続血圧、呼気中二酸化炭素濃度などの実測値を使用するが、データ収集時の被験者の負担が少なく、後に述べる数式モデルの構造を簡略化するための方法として(1)、(2)、(3)の処理をへた生体機能計測データの代表値を使用することも可能である。本実施例では生体機能計測データの代表値として全脳平均値を使用する。なお、病態により正常と異なる信号変化を示す部分が生体機能計測データに含まれている場合には、当該の画素またはチャンネルの計測値を除外して処理を行う。多断面の磁気共鳴画像を用いる場合には断面ごとの撮像タイミングのずれを補正するため、断面毎の撮像タイミングのずれを調整した補間処理を行った後、平均値を算出する。これら測定値に対し、(3)と同様にして直流および超低周波数成分を除去する平均値除去、トレンド除去処理を行う。また、高周波数成分を除去するローパスフィルター処理を行う。なお、データによってはこの処理を省略しても差し支えない。生理学的変動の参照データに遅延がある場合、たとえばサイドストリームとして採取した呼気ガスサンプルを用いて測定した二酸化炭素濃度を使用する場合などには、あらかじめ遅延時間を補正した計測値を入力信号として使用する。生理学的変動の参照データのサンプリングレートが生体機能計測データのそれと異なる場合には両者が同じサンプリングレートとなるようにリサンプリング処理により調節する。
【0012】
(5)数式モデルの生成
利用する数式モデルとしては線形、時不変性モデルに限らず、非線形、時変性モデルを含むあらゆるブラックボックスモデルを使用することができる。生体の活動は一般的には非線形、時変性であるが、安定した生理学的状態、比較的短い時間の範囲では線形、時不変性システムと仮定して取り扱っても大きな誤差は生じない。本実施例ではモデルパラメーターの計算が容易な自己回帰システムモデル(auto−regressive exogenous model; ARX model)を使用して説明する(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著.MATLABによる制御のためのシステム同定,東京電機大学出版局.Ljung L, System Identification Toolbox User’s Guide Version 5, The Mathworks, Inc)。数式を次に示す。
【0013】
【数1】
【0014】
ただし、y(k)はシステムの出力信号、u(k)はシステムの入力信号、nkはむだ時間、w(k)は平均値0、有限な分散を持つ白色雑音、A(q)とB(q)はシフトオペレーターqの多項式である。
画素またはチャンネル毎の計測値を出力、脳機能賦活パラダイムデータと生理学的変動の参照データを入力としてモデルパラメーターの推定を行う。ARXモデルパラメーターの推定には最小二乗推定法を使用するが、これ以外のモデル、たとえば自己回帰移動平均システムモデル(auto−regressive moving average exogenous model; ARMAX model)を使用する場合にはモデルパラメーター推定法として最尤推定法などを使用する(中溝高好著,信号解析とシステム同定,コロナ社)。
【0015】
示適モデル次数はクロスバリデーション、最小実現、情報規範などを用いて選定する(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著,MATLABによる制御のためのシステム同定,東京電機大学出版局)。本実施例では多数の画素またはチャンネルの計測値を処理するにあたって、オペレーターの関与を少なくする方法として情報規範を用いた自動的選定、あるいは自動的に選定されたモデル次数の結果を参考にしたオペレーター選定によった。自動選定に用いる情報規範としては赤池情報量基準などを使用する(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。同様にしてむだ時間nkについても相関法を用いたインパルス応答あるいはクロスバリデーション、情報規範などで選定を行う。なお、データによってはむだ時間を0に固定しても差し支えない。
【0016】
全脳平均値を入力信号として使用した場合には、全脳平均値の変化に先行して信号変化を生ずる脳部位が存在する可能性があるため、全脳平均値の遅延時間の調整を行う。調整量はクロスバリデーションに基づく選定、あるいは既知の通過時間に基づきオペレーターの指示による選定を行う。なお、データによってはこの処理を省略しても差し支えない。
【0017】
局所生理機能の情報抽出のみを目的とする場合には、生理学的変動の参照データを入力とするモデルを使用する。あるいは複数の生理学的変動の参照データを用い、これらの間のフィードバック関係を含めて情報抽出を行う場合は多変量時系列モデルを使用する(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。明らかな入力を特定しえない生理学的変動、たとえば血管運動に起因する信号変化の特徴を抽出する場合は他の生理学的変動の参照データとともに確率的雑音を入力としてモデル化する。あるいは他の生理学的変動の影響が除外できる場合は時系列モデルを使用する。
【0018】
(6)結果の表示
画素、チャンネル毎に推定されたモデルパラメーターは多項式表現、伝達関数表現、零点−極点−ゲイン表現などに形式を変換して表示することができる。あるいは周波数伝達関数として表示することができる。零点、極点などのモデル固有値は複素平面上での分布表示あるいは周波数帯域別の空間的分布を示すマップ表示として示すことができる。
図3はヒトの視覚刺激時の頭部磁気共鳴画像データに本発明を適用し得られたモデルに基づき、脳機能賦活パラダイムに対応するゲインの空間的分布を示すマップである。ゲインの値は画像の輝度により表示している。
図4はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルの固有値(極点)の分布を複素平面上に表示したグラフである。
図5および図6は本発明を適用し得られたモデルの固有値(極点および零点)の周波数帯域別の空間的分布を示すマップである。これらの例においては周波数帯域 0.008〜0.05 Hz、0.05〜0.10 Hz、0.10〜0.15 Hz、0.15〜0.20 Hz、0.20〜0.25 Hzにおける画素毎の固有値の数を左上から右下に向かって配置された画像によって示す。なお最下段右はブランクである。固有値の数は画像の輝度により表示している。図6では、最上段右の周波数帯域 0.05〜0.10 Hzの画像のみが現われている。
確率的雑音の強度は分散あるいは標準偏差、変動係数などに変換して、その分布あるいは空間的分布を示すマップとして表示することができる。
図7はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルの確率的雑音の強度の空間的分布を示すマップである。確率的雑音の強度は画像の輝度により表示している。この例では表示スケールは画素毎の確率的雑音の分散に1を加えた値の常用対数により示している。
モデルに基づいてシミュレーションを行う場合にはインパルス応答、ステップ応答、その他、オペレーターの指定する任意の入力波形に対する応答を表示する。応答は入力開始時を基準とする経時的変化、あるいはオペレーターの指定する経過時間後の信号値の空間的分布を示すマップとして表示することができる。
【0019】
生体機能計測の実測値と入力信号に基づいたシミュレーション結果を比較することで、推定されたモデルの妥当性を検証することができる。
図8はヒトの視覚刺激時の頭部磁気共鳴画像データに本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。この例においては8秒間の強度1のパルス状入力に対する脳内各部位の0.5秒間隔の反応を左上から右下に向かって配置された画像により示す。2段目左端の画像が入力開始時、上から4段目左から4番目の画像が入力終了時に対応する。入力に対応する信号強度の変化量は画像の輝度により表示している。
図9は後頭葉の1画素の信号変化のシミュレーション結果を示すグラフである。横軸は時間経過を入力開始時点を1とする画像番号により表示している。縦軸は信号強度の変化量を示す。参照のため入力波形を重ねて表示している。
図10は上記モデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。この例においては8秒間の強度1のパルス状入力に対する脳内各部位の0.5秒間隔の反応を左上から右下に向かって配置された画像により示す。2段目左端の画像が入力開始時、上から4段目左から4番目の画像が入力終了時に対応する。入力に対応する信号強度の変化量は画像の輝度により表示している。
図11は頭頂葉の1画素の信号変化のシミュレーション結果を示すグラフである。横軸は時間経過を入力開始時点を1とする画像番号により表示している。縦軸は信号強度の変化量を示す。参照のため入力波形を重ねて表示している。
【0020】
本実施例は、脳機能賦活パラダイムとして1種類の入力を使用する例を示したが、複数の入力を用いる脳機能賦活パラダイムに対しても適用することができる。
本実施例は、1断面の磁気共鳴画像を用いているが、多断面のデータでも各断面に対して同様の手順で取り扱うことができる。多断面の場合には断面毎の撮像タイミングのずれを補正して処理を行う。
本実施例では、全データ収集後の一括処理を例にとって説明したが、逐次計算アルゴリズムを使用することにより(中溝高好著,信号解析とシステム同定,コロナ社,足立修一著,MATLABによる制御のためのシステム同定,東京電機大学出版局)、オンラインあるいはリアルタイム処理が可能である。
本実施例ではヒトを対象とした脳機能賦活情報の抽出を例にとって説明したが、本発明は脳以外の臓器、組織およびヒト以外の動物を対象として適用することができる。
【0021】
時系列モデルを使用した場合にはモデルパラメーターの表示または形式を変換した表示、固有値の表示、パワースペクトル、クロススペクトル、コヒーレンスなどの周波数領域の情報に変換した表示、あるいは確率的雑音の強度の表示をすることができる。また、モデルに基づくシミュレーションを表示することができる。多変量自己回帰モデルにおいて、各変量の雑音がお互いに無相関とみなせる場合には、変量間にフィードバックが存在する系であっても、閉ループの周波数応答関数、ノイズ寄与率(またはパワー寄与率)を表示することができる(赤池弘次、中川東一郎共著,ダイナミックシステムの統計的解析と制御[新訂版],サイエンス社.尾崎統、北川源四郎編、時系列解析の方法、朝倉書店)。
図12はヒトの安静時の頭部磁気共鳴画像データと生理学的変動の参照値(脈波、胸郭運動、呼気終末二酸化炭素濃度)に本発明を適用し得られたモデルに基づき、横軸に周波数、縦軸に磁気共鳴信号に対する各変量のノイズ(パワー)寄与率を表示したグラフである。図中、Aは呼気終末二酸化炭素濃度、Bは胸郭運動、Cは脈波、Dは自身(磁気共鳴信号)をそれぞれ示す。
また、本発明を適用し得られたモデルに基づき、周波数帯域を選択して生体機能計測データに対する各変量のノイズ(パワー)寄与率の空間的分布をマップとして表示することができる。
図13は本発明を適用し得られたモデルに基づき、周波数帯域0.98〜1.10 Hzにおける磁気共鳴信号に対する脈波のノイズ(パワー)寄与率の空間的分布を示すマップである。脈波のノイズ(パワー)寄与率は画像の輝度により表示している。
図14は本発明を適用し得られたモデルに基づき、周波数帯域0.24〜0.36 Hzにおける磁気共鳴信号に対する胸郭運動のノイズ(パワー)寄与率の空間的分布を示すマップである。胸郭運動のノイズ(パワー)寄与率は画像の輝度により表示している。
図15は本発明を適用し得られたモデルに基づき、周波数帯域0.04〜0.16 Hzにおける磁気共鳴信号に対する呼気終末二酸化炭素濃度のノイズ(パワー)寄与率の空間的分布を示すマップである。呼気終末二酸化炭素濃度のノイズ(パワー)寄与率は画像の輝度により表示している。
本実施例では変量毎にノイズ(パワー)寄与率の空間的分布を示すマップを作成したが、各変量に色を割り付けることにより1枚のカラーマップとして表示することも可能である。
また、任意の変量に対し入力を与えた場合の応答をシミュレーションとして表示することができる(尾崎統、北川源四郎編、時系列解析の方法、朝倉書店.尾崎統、北川源四郎編、時系列解析の実際I、朝倉書店)。
【0022】
本発明に係る生体機能情報の分別抽出装置としては、例えば、コンピュータを用いて生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出するようにすることができる。分別抽出された生体機能情報はコンピュータに接続された表示装置に表示することができる。
【発明の効果】
本発明では、既知の抽出法である平均値の有意差検定や相関を用いた方法とは異なり、刺激提示あるいは課題遂行に伴う生体機能計測データの信号変化のパターンを仮定する必要がない。また、本発明では統計的検定の前提となる事象の独立性を満足するための平滑化処理を必要としないので生体機能計測データの分解能を劣化させることがない。本発明では独立成分分析とは異なり、分別された信号成分の生理学的由来について推測する必要がない。本発明では、局所神経活動以外の生理学的変動が生体機能計測データに混在する場合でも精度よく局所神経活動の情報を抽出することができる。さらに、局所神経活動以外の種々の生理機能についての情報を抽出することができる。
【図面の簡単な説明】
【図1】生理機能賦活情報の分別抽出手順を示すフローチャートである。
【図2】(a)、(b)は生体機能および生理学的データの収集方法を示す図である。
【図3】本発明を適用し得られたモデルに基づき、脳機能賦活パラダイムに対応するゲインの空間的分布を示すマップである。
【図4】本発明を適用し得られたモデルの固有値(極点)の分布を複素平面上に表示したグラフである。
【図5】本発明を適用し得られたモデルの固有値(極点)の周波数帯域別の空間的分布を示すマップである。
【図6】本発明を適用し得られたモデルの固有値(零点)の周波数帯域別の空間的分布を示すマップである。
【図7】本発明を適用し得られたモデルの確率的雑音の強度の空間的分布を示すマップである。
【図8】本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。
【図9】本発明を適用し得られたモデルに対し、パルス状の入力波形を脳機能賦活パラダイムとして与えたときの1画素の信号変化のシミュレーション結果を示すグラフである。
【図10】本発明を適用し得られたモデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの脳内の各部位の経時的信号変化のシミュレーション結果を示す図である。
【図11】本発明を適用し得られたモデルに対し、パルス状の入力波形を全脳平均値の変動として与えたときの1画素の信号変化のシミュレーション結果を示すグラフである。
【図12】本発明を適用し得られたモデルに基づき、磁気共鳴信号に対する各変量のノイズ(パワー)寄与率を表示したグラフである。
【図13】本発明を適用し得られたモデルに基づき、周波数帯域0.98〜1.10 Hzにおける磁気共鳴信号に対する脈波のノイズ(パワー)寄与率の空間的分布を示すマップである。
【図14】本発明を適用し得られたモデルに基づき、周波数帯域0.24〜0.36 Hzにおける磁気共鳴信号に対する胸郭運動のノイズ(パワー)寄与率の空間的分布を示すマップである。
【図15】本発明を適用し得られたモデルに基づき、周波数帯域0.04〜0.16 Hzにおける磁気共鳴信号に対する呼気終末二酸化炭素濃度のノイズ(パワー)寄与率の空間的分布を示すマップである。
【符号の説明】
1 生体機能データの計測
2 解析領域の選択
3 前処理
4 入力信号の生成
5 数式モデルの生成
6 結果の表示
Claims (13)
- 生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出することを特徴とする生体機能情報の分別抽出方法。
- 抽出される生体機能が中枢神経系の機能であることを特徴とする請求項1記載の生体機能情報の分別抽出方法。
- 抽出される生体機能が局所生理機能であることを特徴とする請求項1記載の生体機能情報の分別抽出方法。
- 数式モデルの構築において、生体機能計測データの代表信号値を使用することを特徴とする請求項1〜3のいずれかに記載の生体機能情報の分別抽出方法。
- 数式モデルの構築において、全身生理機能の測定値を使用することを特徴とする請求項1〜3のいずれかに記載の生体機能情報の分別抽出方法。
- 数式モデルの構築において、確率的雑音を使用することを特徴とする請求項1〜3のいずれかに記載の生体機能情報の分別抽出方法。
- 分別抽出された生体機能情報を表示することを特徴とする請求項1〜6のいずれかに記載の生体機能情報の分別抽出方法。
- 生体機能情報の表示において、モデルに基づくシミュレーションを使用することを特徴とする請求項7記載の生体機能情報の分別抽出方法。
- 生体機能情報の表示において、ノイズ寄与率(パワー寄与率)の分布または空間的情報の表示を使用することを特徴とする請求項7記載の生体機能情報の分別抽出方法。
- 生体機能情報の表示において、モデル固有値の分布または空間的情報の表示を使用することを特徴とする請求項7記載の生体機能情報の分別抽出方法。
- 生体機能情報の表示において、確率的雑音の強度の分布または空間的情報の表示を使用することを特徴とする請求項7記載の生体機能情報の分別抽出方法。
- 生体機能計測データに対して画素またはチャンネル毎に入出力関係を記述する数式モデルを構築し、生体機能計測データの信号に含まれる生理学的起源の異なる構成成分の生体機能情報を分別抽出する手段を備えたことを特徴とする生体機能情報の分別抽出装置。
- 分別抽出された生体機能情報を表示する手段を備えたことを特徴とする請求項12記載の生体機能情報の分別抽出装置。
Priority Applications (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2002265051A JP2004041670A (ja) | 2002-05-17 | 2002-09-11 | 生体機能情報の分別抽出方法およびその装置 |
US10/462,642 US20040049484A1 (en) | 2002-09-11 | 2003-06-17 | Method and apparatus for separating and extracting information on physiological functions |
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2002143021 | 2002-05-17 | ||
JP2002265051A JP2004041670A (ja) | 2002-05-17 | 2002-09-11 | 生体機能情報の分別抽出方法およびその装置 |
Publications (1)
Publication Number | Publication Date |
---|---|
JP2004041670A true JP2004041670A (ja) | 2004-02-12 |
Family
ID=31719475
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2002265051A Pending JP2004041670A (ja) | 2002-05-17 | 2002-09-11 | 生体機能情報の分別抽出方法およびその装置 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP2004041670A (ja) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2009273539A (ja) * | 2008-05-13 | 2009-11-26 | Research Organization Of Information & Systems | イメージングデータからの状態変化抽出方法、状態変化の可視化システム、および、コンピュータプログラム |
JP2014507208A (ja) * | 2011-01-19 | 2014-03-27 | ドルフィン テクノロジーズ オサケ ユキチュア | 心臓血管脈拍波動の視覚化の方法および装置 |
CN108366752A (zh) * | 2015-11-24 | 2018-08-03 | 株式会社国际电气通信基础技术研究所 | 脑活动分析装置、脑活动分析方法、程序以及生物标记物装置 |
WO2023021671A1 (ja) * | 2021-08-19 | 2023-02-23 | 国立大学法人大阪大学 | 壁厚み推定方法、壁厚み推定装置及び壁厚み推定システム |
-
2002
- 2002-09-11 JP JP2002265051A patent/JP2004041670A/ja active Pending
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2009273539A (ja) * | 2008-05-13 | 2009-11-26 | Research Organization Of Information & Systems | イメージングデータからの状態変化抽出方法、状態変化の可視化システム、および、コンピュータプログラム |
JP2014507208A (ja) * | 2011-01-19 | 2014-03-27 | ドルフィン テクノロジーズ オサケ ユキチュア | 心臓血管脈拍波動の視覚化の方法および装置 |
CN108366752A (zh) * | 2015-11-24 | 2018-08-03 | 株式会社国际电气通信基础技术研究所 | 脑活动分析装置、脑活动分析方法、程序以及生物标记物装置 |
WO2023021671A1 (ja) * | 2021-08-19 | 2023-02-23 | 国立大学法人大阪大学 | 壁厚み推定方法、壁厚み推定装置及び壁厚み推定システム |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Agrawal et al. | Model-based physiological noise removal in fast fMRI | |
Kortelainen et al. | Sleep staging based on signals acquired through bed sensor | |
Kassinopoulos et al. | Identification of physiological response functions to correct for fluctuations in resting-state fMRI related to heart rate and respiration | |
Mayeli et al. | Real-time EEG artifact correction during fMRI using ICA | |
Mainardi et al. | Pole-tracking algorithms for the extraction of time-variant heart rate variability spectral parameters | |
DE102017112819A1 (de) | System und Verfahren zum Vorsehen einer Echtzeitsignal-Segmentierung und Kontrollpunktausrichtungsrahmen | |
JP7310048B2 (ja) | 診断支援プログラム | |
Falahpour et al. | Subject specific BOLD fMRI respiratory and cardiac response functions obtained from global signal | |
US10898143B2 (en) | Quantitative mapping of cerebrovascular reactivity using resting-state functional magnetic resonance imaging | |
EP4193912A1 (en) | Method and system for use in monitoring neural activity in a subject's brain | |
US10646165B2 (en) | Removing eletrophysicologic artifacts from a magnetic resonance imaging system | |
Kar et al. | Effect of sleep deprivation on functional connectivity of EEG channels | |
US10349896B2 (en) | Epsilon-tube filter for blunt noise removal | |
US20220183561A1 (en) | Generating imaging-based neurological state biomarkers and estimating cerebrospinal fluid (csf) dynamics based on coupled neural and csf oscillations during sleep | |
Nasrolahzadeh et al. | Analysis of heart rate signals during meditation using visibility graph complexity | |
Prokopiou et al. | Modeling the hemodynamic response function using EEG-fMRI data during eyes-open resting-state conditions and motor task execution | |
US20040049484A1 (en) | Method and apparatus for separating and extracting information on physiological functions | |
JP2004041670A (ja) | 生体機能情報の分別抽出方法およびその装置 | |
KR101919907B1 (ko) | 다중 신경생리신호 기반 사용자 간 상호작용 모니터링 장치 및 방법 | |
de la Cruz et al. | Impact of the heart rate on the shape of the cardiac response function | |
Kurbako et al. | Mathematical models of the electrocardiogram and photoplethysmogram signals to test methods for detection of synchronization between physiological oscillatory processes | |
Keenan et al. | Adaptive filtering of heart rate signals for an improved measure of cardiac autonomic control | |
Khouaja et al. | Enhancing eeg surface resolution by using a combination of kalman filter and interpolation method | |
JP2814969B2 (ja) | 活動部位検出装置 | |
Borovkova et al. | Experimental observation of Arnold tongues in the analysis of the signal from contour of the autonomous regulation of heart rate and respiration |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20050815 |
|
A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20071109 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20071120 |
|
A02 | Decision of refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A02 Effective date: 20080318 |