(57)【要約】【課題】  短い解析時間で精度の高い非散乱体を含む不均一散乱体の成分濃度等の情報を同定する。【解決手段】  光パルス応答を用いて不均一散乱体に含まれる成分の情報を同定する成分情報同定装置１０は、不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定部１１と、散乱体部分に対し複数のセルを設定する離散化部１２と、光拡散方程式の微分成分を差分成分として置換え、各セルにおける光伝搬の解析結果を導出する解析結果導出部１３と、非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として非散乱体部分を伝搬する光を上記光拡散方程式の等価波源として組み込む等価波源組込部１４と、光伝搬の解析結果から光パルス応答の理論値を導出する理論値導出部１５と、光パルス応答の理論値と測定値とから光吸収係数を同定する光吸収係数同定部１６とを備える。【選択図】  図３
2005221488.tif 000002
【特許請求の範囲】
【請求項１】
不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報を同定する成分情報同定方法において、
不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、
前記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、
光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、前記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、
前記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、前記光拡散方程式の等価波源として組み込む等価波源組込ステップと、
前記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、
前記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップと
を備える成分情報同定方法。
【請求項２】
前記光吸収係数同定ステップにおいて同定された光吸収係数から、前記成分の濃度を同定する成分濃度同定ステップを
更に備える請求項１に記載の成分情報同定方法。
【請求項３】
前記光吸収係数同定ステップにおいて、前記光吸収係数を同定すると共に、前記比較の結果を参照することにより光散乱係数を同定する
ことを特徴とする請求項１に記載の成分情報同定方法。
【請求項４】
前記離散化ステップにおいて、前記散乱体部分に対し格子状に分割されるような前記複数のセルを設定し、
前記解析結果導出ステップにおいて、前記光拡散方程式による導出を、前記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う
ことを特徴とする請求項１〜３の何れか一項に記載の成分情報同定方法。
【請求項５】
前記光吸収係数同定ステップにおいて、前記理論値と前記実測値との自乗残差最小誤差を最小とする基準に基づき、前記光吸収係数を同定する
ことを特徴とする請求項１〜４の何れか一項に記載の成分情報同定方法。
【請求項６】
不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報を同定する成分情報同定装置において、
不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定手段と、
前記設定手段により設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化手段と、
光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、前記離散化手段により設定された各セルにおける光伝搬の解析結果を導出する解析結果導出手段と、
前記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、前記光拡散方程式の等価波源として組み込む等価波源組込手段と、
前記解析結果導出手段により導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出手段と、
前記理論値導出手段により導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定手段と
を備える成分情報同定装置。
【請求項７】
前記光吸収係数同定手段により同定された光吸収係数から、前記成分の濃度を同定する成分濃度同定手段を
更に備える請求項６に記載の成分情報同定装置。
【請求項８】
前記光吸収係数同定手段は、前記光吸収係数を同定すると共に、前記比較の結果を参照することにより光散乱係数を同定する
ことを特徴とする請求項６に記載の成分情報同定装置。
【請求項９】
前記離散化手段は、前記散乱体部分に対し格子状に分割されるような前記複数のセルを設定し、
前記解析結果導出手段は、前記光拡散方程式による導出を、前記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う
ことを特徴とする請求項６〜８の何れか一項に記載の成分情報同定装置。
【請求項１０】
前記光吸収係数同定手段は、前記理論値と前記実測値との自乗残差最小誤差を最小とする基準に基づき、前記光吸収係数を同定する
ことを特徴とする請求項６〜９の何れか一項に記載の成分情報同定装置。
【請求項１１】
前記不均一散乱体の光パルス応答を実測する実測手段を更に備え、
前記光吸収係数同定手段は、前記理論値と、前記実測手段により実測された光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する
ことを特徴とする請求項６〜１０の何れか一項に記載の成分情報同定装置。
【請求項１２】
不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報をコンピュータに同定させる成分情報同定プログラムにおいて、
不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、
前記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、
光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、前記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、
前記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、前記光拡散方程式の等価波源として組み込む等価波源組込ステップと、
前記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、
前記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップと
を前記コンピュータに実行させる成分情報同定プログラム。
【請求項１３】
前記光吸収係数同定ステップにおいて同定された光吸収係数から、前記成分の濃度を同定する成分濃度同定ステップを
前記コンピュータに更に実行させる請求項１２に記載の成分情報同定プログラム。
【請求項１４】
前記光吸収係数同定ステップにおいて、前記光吸収係数を同定すると共に、前記比較の結果を参照することにより光散乱係数を同定する
ことを特徴とする請求項１２に記載の成分情報同定プログラム。
【請求項１５】
前記離散化ステップにおいて、前記散乱体部分に対し格子状に分割されるような前記複数のセルを設定し、
前記解析結果導出ステップにおいて、前記光拡散方程式による導出を、前記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う
ことを特徴とする請求項１２〜１４の何れか一項に記載の成分情報同定プログラム。
【請求項１６】
前記光吸収係数同定ステップにおいて、前記理論値と前記実測値との自乗残差最小誤差を最小とする基準に基づき、前記光吸収係数を同定する
ことを特徴とする請求項１２〜１５の何れか一項に記載の成分情報同定プログラム。
【請求項１７】
不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報をコンピュータに同定させる成分情報同定プログラムを記録したコンピュータ読取可能な記録媒体において、
不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、
前記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、
光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、前記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、
前記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、前記光拡散方程式の等価波源として組み込む等価波源組込ステップと、
前記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、
前記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップと
を前記コンピュータに実行させる成分情報同定プログラムを記録したコンピュータ読取可能な記録媒体。
【請求項１８】
前記光吸収係数同定ステップにおいて同定された光吸収係数から、前記成分の濃度を同定する成分濃度同定ステップを
前記コンピュータに更に実行させる請求項１７に記載の記録媒体。
【請求項１９】
前記光吸収係数同定ステップにおいて、前記光吸収係数を同定すると共に、前記比較の結果を参照することにより光散乱係数を同定する
ことを特徴とする請求項１７に記載の記録媒体。
【請求項２０】
前記離散化ステップにおいて、前記散乱体部分に対し格子状に分割されるような前記複数のセルを設定し、
前記解析結果導出ステップにおいて、前記光拡散方程式による導出を、前記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う
ことを特徴とする請求項１７〜１９の何れか一項に記載の記録媒体。
【請求項２１】
前記光吸収係数同定ステップにおいて、前記理論値と前記実測値との自乗残差最小誤差を最小とする基準に基づき、前記光吸収係数を同定する
ことを特徴とする請求項１７〜２０の何れか一項に記載の記録媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報を同定する成分情報同定方法、成分情報同定装置、当該情報をコンピュータに同定させる成分情報同定プログラム、及び当該成分情報同定プログラムを記憶した記憶媒体に関する。
【背景技術】
【０００２】
光トモグラフィにより生体中の成分濃度を同定する手法として、生体に光を入射したときの光応答を測定し、その測定値から光学パラメータを推定し、その光学パラメータを用いて同定する方法がある。例えば特許文献１では、光源と検出器を対称に配置して、その差信号をとることにより生体の局部的な吸収係数変化を推定することにより、生体の血中成分濃度を同定する技術が言及されている。
【特許文献１】特開２００２−３６９８１３号公報
【発明の開示】
【発明が解決しようとする課題】
【０００３】
しかしながら、上記の方法では連続光を用いているため、同定精度及び空間分解能等に改善の余地がある。連続光でなくパルス光を用いる方法もあり、パルス光を用いる方法は成分濃度の絶対値の推定が可能であり光吸収係数及び光散乱係数を個別に推定できるので、連続光を用いる方法より同定精度及び空間分解能等が優れている。しかし、パルス光を用いる方法は、光パルス応答の演算等の点で解析が困難である。
【０００４】
また、成分濃度を同定する対象として、人体の頭部の血中に含まれる酸化ヘモグロビン及び還元ヘモグロビンが挙げられる。人体の頭部は、散乱体である頭蓋骨と灰／白質との間に非散乱体である脳髄液が介在する不均一散乱体である。上記パルス光を用いる解析では、通常、散乱体のみで構成されている生体部位を想定したものであり、上記人体の頭部のような非散乱体を含んだ不均一散乱体を精度よく解析することができない。非散乱体を含む不均一散乱体を解析する方法も提案されているが、何れの方法も、解析時間や精度の点で十分なものであるとはいえない。
【０００５】
本発明は、上記課題を解決するために成されたものであり、短い解析時間で精度の高い非散乱体を含む不均一散乱体の成分濃度等の情報を同定することができる成分情報同定方法、成分情報同定装置、成分情報同定プログラム、及び成分情報同定プログラムを記録したコンピュータ読取可能な記録媒体を提供することを目的とする。
【課題を解決するための手段】
【０００６】
上記目的を達成するために、本発明に係る成分情報同定方法は、不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報を同定する成分情報同定方法において、不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、上記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、上記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、上記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、上記光拡散方程式の等価波源として組み込む等価波源組込ステップと、上記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、上記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップとを備えることを特徴とする。
【０００７】
ここで、不均一散乱体に含まれる成分の情報には、例えば、光吸収係数等の光学パラメータや成分濃度等が相当する。
【０００８】
本発明に係る成分情報同定方法では、設定ステップにおいて、解析対象である不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する。以降は、設定された散乱体部分及び非散乱体部分毎に解析を行う。
【０００９】
続いて、離散化ステップにおいて、上記設定された設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する。
【００１０】
続いて、解析結果ステップにおいて、光拡散方程式を用い演算することにより、上記設定された各セルにおける光伝搬の解析結果を導出する。
【００１１】
一方、等価波源組込ステップにおいて、非散乱体部分を伝搬する光の光強度を算出し、光拡散方程式の等価波源として組み込む。
【００１２】
上記解析結果ステップでの各セルにおける光伝搬の解析結果の導出は、上記等価波源組込ステップでの上記光伝搬の解析結果を利用した光強度の算出及び等価波源の組込行いながら、セルごとに順次行われる。
【００１３】
各セルにおける光伝搬の解析結果が導出された後、理論値導出ステップにおいて、上記解析結果を演算することにより、光パルス応答の理論値を導出する。
【００１４】
続いて、光吸収係数同定ステップにおいて、上記理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する。
【００１５】
上記により、本発明に係る成分情報同定方法によれば、上記散乱体を離散化して単純な代数計算のみで解析することにより解析時間を短くすることができ、非散乱体を含む不均一散乱体の光パルス応答の推定を可能としていることにより精度の高い非散乱体を含む不均一散乱体に含まれる成分の情報を同定することができる。
【００１６】
また、本発明に係る成分情報同定方法は、上記光吸収係数同定ステップにおいて同定された光吸収係数から、上記成分の濃度を同定する成分濃度同定ステップを更に備える構成とすることが望ましい。
【００１７】
上記により、本発明に係る成分情報同定方法は、成分の濃度を同定することができる。
【００１８】
また、本発明に係る成分情報同定方法は、上記光吸収係数同定ステップにおいて、上記光吸収係数を同定すると共に、上記比較の結果を参照することにより光散乱係数を同定する構成とすることが望ましい。
【００１９】
上記により、本発明に係る成分情報同定方法は、光散乱係数を同定することができる。
【００２０】
また、本発明に係る成分情報同定方法は、より具体的には、上記離散化ステップにおいて、上記散乱体部分に対し格子状に分割されるような上記複数のセルを設定し、上記解析結果導出ステップにおいて、上記光拡散方程式による導出を、上記光拡散方程式を時間的に離散化し、ＦＤＴＤ（Finite Difference Time Domain）法により行う構成とすることが望ましい。
【００２１】
また、本発明に係る成分情報同定方法は、より具体的には、上記光吸収係数同定ステップにおいて、上記理論値と上記実測値との自乗残差最小誤差を最小とする基準に基づき、上記光吸収係数を同定する構成とすることが望ましい。
【００２２】
本発明は、以下に述べるような成分情報同定装置、プログラム、及び記憶媒体に係る発明として捉えることもでき、同様の作用及び効果を奏する。
【００２３】
本発明に係る成分情報同定装置は、不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報を同定する成分情報同定装置において、不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定手段と、上記設定手段により設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化手段と、光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、上記離散化手段により設定された各セルにおける光伝搬の解析結果を導出する解析結果導出手段と、上記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、上記光拡散方程式の等価波源として組み込む等価波源組込手段と、上記解析結果導出手段により導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出手段と、上記理論値導出手段により導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定手段とを備えることを特徴とする。
【００２４】
また、本発明に係る成分情報同定装置は、上記光吸収係数同定手段により同定された光吸収係数から、上記成分の濃度を同定する成分濃度同定手段を更に備える構成とすることが望ましい。
【００２５】
また、本発明に係る成分情報同定装置において、上記光吸収係数同定手段は、上記光吸収係数を同定すると共に、上記比較の結果を参照することにより光散乱係数を同定する構成とすることが望ましい。
【００２６】
また、本発明に係る成分情報同定装置において、より具体的には、上記離散化手段は、上記散乱体部分に対し格子状に分割されるような上記複数のセルを設定し、上記解析結果導出手段は、上記光拡散方程式による導出を、上記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う構成とすることが望ましい。
【００２７】
また、本発明に係る成分情報同定装置において、より具体的には、上記光吸収係数同定手段は、上記理論値と上記実測値との自乗残差最小誤差を最小とする基準に基づき、上記光吸収係数を同定する構成とすることが望ましい。
【００２８】
また、本発明に係る成分情報同定装置は、上記不均一散乱体の光パルス応答を実測する実測手段を更に備え、上記光吸収係数同定手段は、上記理論値と、上記実測手段により実測された光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する構成とすることが望ましい。
【００２９】
本発明に係る成分情報同定プログラムは、不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報をコンピュータに同定させる成分情報同定プログラムにおいて、不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、上記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、上記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、上記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、上記光拡散方程式の等価波源として組み込む等価波源組込ステップと、上記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、上記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップとを上記コンピュータに実行させることを特徴とする。
【００３０】
また、本発明に係る成分情報同定プログラムは、上記光吸収係数同定ステップにおいて同定された光吸収係数から、上記成分の濃度を同定する成分濃度同定ステップを上記コンピュータに更に実行させる構成とすることが望ましい。
【００３１】
また、本発明に係る成分情報同定プログラムは、上記光吸収係数同定ステップにおいて、上記光吸収係数を同定すると共に、上記比較の結果を参照することにより光散乱係数を同定する構成とすることが望ましい。
【００３２】
また、本発明に係る成分情報同定プログラムは、より具体的には、上記離散化ステップにおいて、上記散乱体部分に対し格子状に分割されるような上記複数のセルを設定し、上記解析結果導出ステップにおいて、上記光拡散方程式による導出を、上記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う構成とすることが望ましい。
【００３３】
また、本発明に係る成分情報同定プログラムは、より具体的には、上記光吸収係数同定ステップにおいて、上記理論値と上記実測値との自乗残差最小誤差を最小とする基準に基づき、上記光吸収係数を同定する構成とすることが望ましい。
【００３４】
本発明に係る記録媒体は、不均一散乱体に所定の入射位置からパルス光を入射したとき出射位置での光パルス応答を用いて、当該不均一散乱体に含まれる成分の情報をコンピュータに同定させる成分情報同定プログラムを記録したコンピュータ読取可能な記録媒体において、不均一散乱体を散乱体部分及び非散乱体部分の何れかに分けて設定する設定ステップと、上記設定ステップにおいて設定された散乱体部分に対し複数のセルを設定し、空間的に離散化する離散化ステップと、光拡散方程式を用い、当該光拡散方程式の微分成分を差分成分として置換え、演算することにより、上記離散化ステップにおいて設定された各セルにおける光伝搬の解析結果を導出する解析結果導出ステップと、上記非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変として、当該非散乱体部分を伝搬する光の光強度を算出し、上記光拡散方程式の等価波源として組み込む等価波源組込ステップと、上記解析結果導出ステップにおいて導出された各セルにおける光伝搬の解析結果を演算することにより、光パルス応答の理論値を導出する理論値導出ステップと、上記理論値導出ステップにおいて導出された光パルス応答の理論値と、光パルス応答の実測値とを比較し、当該比較の結果を参照することにより、光吸収係数を同定する光吸収係数同定ステップとを上記コンピュータに実行させる成分情報同定プログラムを記録したことを特徴とする。
【００３５】
また、本発明に係る記録媒体は、上記光吸収係数同定ステップにおいて同定された光吸収係数から、上記成分の濃度を同定する成分濃度同定ステップを上記コンピュータに更に実行させる構成とすることが望ましい。
【００３６】
また、本発明に係る記録媒体は、上記光吸収係数同定ステップにおいて、上記光吸収係数を同定すると共に、上記比較の結果を参照することにより光散乱係数を同定する構成とすることが望ましい。
【００３７】
また、本発明に係る記録媒体は、より具体的には、上記離散化ステップにおいて、上記散乱体部分に対し格子状に分割されるような上記複数のセルを設定し、上記解析結果導出ステップにおいて、上記光拡散方程式による導出を、上記光拡散方程式を時間的に離散化し、ＦＤＴＤ法により行う構成とすることが望ましい。
【００３８】
また、本発明に係る記録媒体は、より具体的には、上記光吸収係数同定ステップにおいて、上記理論値と上記実測値との自乗残差最小誤差を最小とする基準に基づき、上記光吸収係数を同定する構成とすることが望ましい。
【発明の効果】
【００３９】
上記により、本発明によれば、解析時間を短くすることができ、精度の高い非散乱体を含む不均一散乱体に含まれる成分の情報を同定することができる。
【発明を実施するための最良の形態】
【００４０】
以下、図面を参照しながら本発明の実施形態について説明する。
【００４１】
図１は、本発明に係る成分情報同定装置を用いた成分情報を同定するシステムの構成を示す。本実施形態では、人頭部の血中酸素濃度分布の測定に適用するものとする。血中酸素濃度として具体的には酸化ヘモグロビン及び還元ヘモグロビン濃度等の情報を測定する。
【００４２】
図１に示すように、本システムは成分情報同定装置１０と、光パルス測定装置２０とを含んで構成され、その測定対象は人頭部３０である。人頭部３０は、頭蓋骨３１、脳髄液３２及び灰／白質３３を含んで構成されている。頭蓋骨３１及び灰／白質３３は、散乱体であり、脳髄液３２は非散乱体であるので、人頭部３０は非散乱体を含む不均一散乱体である。なお、測定対象は上記人頭部３０に限定されるものではなく、不均一散乱体であればよい。
【００４３】
光パルス測定装置２０は、人頭部３０にパルス光を入射し、出射位置での光パルス応答を測定する装置である。図１に示すように、光パルス測定装置２０は、パルス発生器２１と、レーザ光源２２ａ，２２ｂと、光スイッチ２３と、光検出器２４と、波形観測装置２５と、光ガイド２６ａ，２６ｂとを備える。
【００４４】
パルス発生器２１は、レーザ光源２２ａ，２２ｂを駆動させる。レーザ光源２２ａ，２２ｂは、それぞれ互いに異なる波長のパルス光であるレーザ光を発生させる。当該パルス光の波長は、図２のグラフに示す酸化ヘモグロビンＨｂ及び還元ヘモグロビンＨｂＯ２の光吸収係数の波長依存性を考慮して、例えばλ１＝７５０ｎｍ、λ２＝８５０ｎｍとする。光スイッチ２３は、上記２種類のパルス光の波長を選択する。光スイッチ２３により波長を選択されたパルス光は光ガイド２６ａを通り人頭部３０に入射する。人頭部３０内で散乱され、出射位置に到達した光は光ガイド２６ｂを通り、光検出器２４により電気信号に変換される。光検出器２４には、光電子増倍管又はＡＰＤ（Avalanche Photodiode）等を用いる。上記電気信号に変換された光は、波形観測装置２５により、波形を観測される。波形観測装置２５には、オシロスコープ、ストリークカメラ又は時間相関フォトンカウンティング装置等を用いる。観測された波形データは、光パルス応答のデータとして成分情報同定装置１０に送信される。
【００４５】
成分情報同定装置１０は、光パルス測定装置２０から光パルス応答のデータを受信し、人頭部３０の血中酸素濃度分布に関する各種の情報を同定する。図３に示すように、成分情報同定装置１０は、設定部１１と、離散化部１２と、解析結果導出部１３と、等価波源組込部１４と、理論値導出部１５と、光吸収係数同定部１６と、成分濃度同定部１７とを備える。なお、各構成要素についての機能は、処理の説明と共に後述する。
【００４６】
図４に、成分情報同定装置１０のハードウェア構成を示す。上記の各構成要素によるソフトウェア的処理機能、及び成分情報同定装置１０によって行われる成分情報同定方法における各ステップに対応する処理は、これらのハードウェア構成によって行われる。
【００４７】
ＣＰＵ１０１は、ソフトウェアプログラムを実行することにより、各装置の制御やデータの計算、加工を行い、後述する成分情報同定装置１０の各構成要素の機能を司る。ＣＰＵ１０１には、成分情報同定装置１０の処理動作に必要な各ソフトウェアプログラム等が記憶されているＲＯＭ１０２と、プログラム実行中に一時的にデータが記憶されるＲＡＭ１０３とが接続されている。また、このＣＰＵ１０１には、ハードディスク等の外部記憶装置１０４が接続されている。
【００４８】
これらのＣＰＵ１０１等に対して、使用者が適宜操作及び参照が可能なように、入力装置１０５と、表示装置１０６とが接続されて成分情報同定装置１０が構成されている。入力装置１０５としては、例えばマウス等のポインティングデバイスやキーボード等が用いられる。また、表示装置１０６としては、例えばＣＲＴディスプレイや液晶ディスプレイなどが用いられる。また、他の装置、例えば光パルス測定装置２０、との間でデータの送受信をするため、通信インターフェース１０７が設けられている。
【００４９】
なお、成分情報同定装置１０と、光パルス測定装置２０とは別構成としているが、光パルス測定装置２０の各構成要素が成分情報同定装置１０に含まれている構成としてもよい。また、成分情報同定装置１０は、光パルス測定装置２０から直接データを送信される形態でなく、例えばフレキシブルディスク等の媒体を介したデータの取得等、何らかの手段により別途、データを取得する形態としてもよい。
【００５０】
以下、図４のフローチャートを用いて、本実施形態に係る成分情報同定装置１０において実行される処理を説明する。
【００５１】
まず、設定部１１が、測定対象である人頭部３０を散乱体部分及び非散乱体部分に分けて設定する（Ｓ０１、設定ステップ）。即ち、人頭部３０の頭蓋骨３１及び灰／白質３３の部分を散乱体部分と設定し、人頭部３０の脳髄液３２の部分を非散乱体部分と設定する。これは、散乱体部分と非散乱体部分とを別々の手法で分析するためである。
【００５２】
続いて、離散化部１２が、上記設定された散乱体部分に対し、格子状（Ｙｅｅ格子）に分割されるような複数のセルを設定し、空間的に離散化する（Ｓ０２、離散化ステップ）。これは、以降の散乱体部分の分析を各セル単位で行うためである。上記各セルは、本実施形態では立方体とする。本実施形態においては、その一辺は０．５ｍｍ程度とする。この程度のサイズで十分な精度が得ることができる。また、各セル、並びに上記光パルス測定装置２０で測定される波長λ１及びλ２における光学パラメータ、即ち光吸収係数μａ及び等価散乱係数μｓ´の初期値を、文献等に報告されている人頭部３０の各セルの平均的な値を用いて設定する（Ｓ０３）。なお、光散乱係数μｓは、非等方散乱パラメータをｇとしたときに、等価散乱係数μｓ´とμｓ´＝（１−ｇ）μｓの関係である。
【００５３】
以下、解析結果導出部１３及び等価波源組込部１４による光伝搬の解析結果の導出について説明する。解析結果導出部１３は、光拡散方程式を用いて、散乱体部分の分割されたセル及び時間間隔ごとに順次、演算を行う。等価波源組込部１４は、上記セルごとの演算結果を用い、非散乱体部分における光強度の演算を行い、当該演算結果を上記光拡散方程式に組み入れる。即ち、解析結果導出部１３による各セル及び時間間隔毎の演算、並びに等価波源組込部１４による非散乱体部分での光強度の拡散方程式への組込を順次行っていくことにより、光伝搬の解析結果を導出する（Ｓ０４）。
【００５４】
まず、解析結果導出部１３による光拡散方程式を用いた、上記各セルにおける光伝搬の解析結果の導出について説明する。本実施形態における上記光伝搬の解析結果とは、各セルにおける拡散放射流速、放射発散度、及びこれらの値から求められる光強度等をいう。ここで、光拡散方程式は、一般には以下の２つの式で与えられ、当該光拡散方程式を解くことにより散乱体の光伝搬の解析結果が導出される（Ｓ０４ａ、解析結果導出ステップ）。
【数１】
2005221488.tif 000003
【数２】
2005221488.tif 000004
上式において、φ（ｒ，ｔ）は放射発散度で、位置ｒ及び時間ｔの関数である。また、Ｊ（ｒ，ｓ，ｔ）は拡散放射流速で、位置ｒ、方向を示す単位ベクトルｓ、及び時間ｔの関数である。である。また、輸送減衰定数μｔｒは、μｔｒ＝μａ＋（１−ｇ）μｓにより定義される。また、全減衰係数μｔはμｔ＝μａ＋μｓにより定義される。また、ｃは散乱体中の光速である。また、等価波源関数εｒｉ（ｒ，ｓ，ｔ）は、入射光強度Ｉｒｉ（ｒ，ｓ，ｔ）と散乱位相関数ｐ（ｓ，ｓ´）との積を全立体角で積分する次式で与えられる。なお、本文中でのｓは、ハットｓを表すものとする。
【数３】
2005221488.tif 000005
【００５５】
上記光拡散方程式の放射発散度及び拡散放射流速をセルの頂点及び中点に配置し、微分を差分で置き換えると、以下のＦＤＴＤ解が得られる。
【数４】
2005221488.tif 000006
【数５】
2005221488.tif 000007
【数６】
2005221488.tif 000008
【数７】
2005221488.tif 000009
上式において、ｘ，ｙ，ｚは空間の座標軸を示し、ｚはパルス光の入射方向の軸であるとする。また、（ｉ，ｊ，ｋ）はセルを形成している格子座標である。また、Ｆ０（ρ）は、入射光の（ｘ，ｙ）面における光強度分布である。ｆ（ｔ）は、入射光の時間波形である。また、Ｊｘｎ＋１／２（ｉ＋１／２，ｊ，ｋ）は時刻ｔ＝（ｎ＋１／２）Δｔにおける格子座標（ｉ＋１／２，ｊ，ｋ）での拡散放射流速のｘ成分である。また、Ｊｙｎ＋１／２（ｉ，ｊ＋１／２，ｋ）は時刻ｔ＝（ｎ＋１／２）Δｔにおける格子座標（ｉ，ｊ＋１／２，ｋ）での拡散放射流速のｙ成分である。また、Ｊｚｎ＋１／２（ｉ，ｊ，ｋ＋１／２）は時刻ｔ＝（ｎ＋１／２）Δｔにおける格子座標（ｉ，ｊ，ｋ＋１／２）での拡散放射流速のｚ成分である。φｎ（ｉ，ｊ，ｋ）は、時刻ｔ＝（ｎ＋１／２）Δｔにおける格子座標（ｉ，ｊ，ｋ）での放射発散度である。なお、上記においてΔｔは時間間隔を示す。また、ｚ０は散乱体に入射したフォトンが最初に散乱される深さであり、ｚ０＝１／μｓ´で与えられる。また、上記式中の係数は、Ｓ０２において設定された各セルにおける光学パラメータ、各セルのサイズΔｚ、及び時間間隔Δｔにより、各格子座標ｒにおいて、次式のように定義される。
【数８】
2005221488.tif 000010
なお、本実施形態において時間間隔Δｔは、後述する式（２４）からピコ秒オーダの値とするのがよい。
【００５６】
本実施形態において、解析結果導出部１３は微分成分が差分成分で置換えられた上記の光拡散方程式を用いる。即ち解析結果導出部１３は、φ（ｒ，ｔ）をｎΔｔの時間ごとに、Ｊ（ｒ，ｓ，ｔ）を（ｎ＋１／２）Δｔの時間ごとに、順次計算する（Ｓ０４ａ、解析結果導出ステップ）。上式から分かるように、上記の順次計算は、単純な代数計算により、行うことができる。
【００５７】
なお、上記の光拡散方程式（４）〜（８）は、パルス光の入射方向をｚ軸方向として、微分成分を差分成分として置換えたものである。もし、パルス光の入射方向を座標空間上で上記と異なる方向とすれば、その座標空間上でのパルス光の入射方向に応じた式となる。
【００５８】
ここで、散乱体表面ｒｐでの光パルス応答Ｐ（ｒｐ，ｎΔｔ）は次式で与えられる。
【数９】
2005221488.tif 000011
【００５９】
上記説明した光拡散方程式は、非散乱体部分での光伝搬が考慮されていない。以下に、非散乱体部分での光伝搬についての考え方と、等価波源組込部１４による等価波源の計算と、計算された等価波源の上記光拡散方程式への組み込みを説明する。
【００６０】
非散乱体を含む不均一散乱体の光パルス応答は、非散乱体中では光散乱方程式は成立しないので、図６に示すように散乱体間の光結合を考慮することにより推定可能である。任意の時刻ｔにおける非散乱体と散乱体との境界ｒｂにおける光強度Ｉ（ｒｂ，ｓθ，ｔ）は、光拡散近似を用いて次式で与えられる。
【数１０】
2005221488.tif 000012
上式において、ｓθは角度を示す単位ベクトルである。非散乱体部分での光伝搬が自由空間伝播であり光強度空間不変とすると、上記の光強度Ｉ（ｒｂ，ｓθ，ｔ）を有する非散乱光による、非散乱体を挟んだ、非散乱体と散乱体との対面の境界ｒｓにおいて励起される散乱光を以下の式により等価波源ε（ｒｓ，ｓ，ｔ）として導出することができる。
【数１１】
2005221488.tif 000013
【数１２】
2005221488.tif 000014
上式において、ｒ＝｜ｒｓ−ｒｂ｜である。また、θは（ｒｓ−ｒｂ）とｚ軸のなす角、ΔΩは散乱体中のｒｓからみた格子座標ｒｂに対応するセルの立体角であり、セル内の光強度は等しいと仮定している。また、αｃは非散乱体の光吸収係数、ｌｃは非散乱体中でのＩ（ｒｂ，ｓθ，ｔ）の光路長である。ｔ０はｒｂからｒｓへの伝搬時間である。非散乱体中を伝搬した非散乱光は境界面からｚ０以上離れたセルで最初に散乱光を励起すると仮定した。即ち、散乱体入射後、最初の距離ｌｓ１（＝ｚ０）では吸収のみが生じ、残りの距離ｌｓ２では散乱により減衰すると仮定した。上記の式（１１）及び式（１２）で与えられる等価波源による拡散放射流速各成分及び放射発散度の変化は、上記境界全てのｒｂの波源を加算して、次の２式で与えられる。なお、式（１１）、式（１２）及び図６では簡単のため、非散乱体及び散乱体の屈折率が等しい場合の計算プロセスを示している。両者の屈折率が異なる場合は、境界でスネルの法則を適用することにより、容易に等価波源を導出可能である。この場合、境界での反射光による等価波源も考慮する必要がある。
【数１３】
2005221488.tif 000015
【数１４】
2005221488.tif 000016
式（１３）及び式（１４）において、Ｂは非散乱体と散乱体との境界面のセル頂点の集合である。図６では、頭蓋骨３１からの非散乱光のみを示しているが、灰／白質３３からの非散乱光による等価波源も同様に考慮する。等価波源組込部１４は、これらの等価波源を計算し、光拡散方程式に組み込む（Ｓ０４ｂ、等価波源組込ステップ）。当該組み込みがなされた上記光拡散方程式を用いて、解析結果導出部１３は、散乱体中の光伝搬の解析結果を導出する（Ｓ０４ａ、解析結果導出ステップ）。なお、式（１２）においてμｔが他の数に比べ相対的に大きな値をもつので、等価波源は散乱体部分と非散乱体部分との境界から１０階層程度のセルで組み入れられれば、実際の計算には十分である。
【００６１】
また、式（４）〜（７）とその係数を与える式（８）、及び式（１３），（１４）で与えられる光結合項を加え、各セルの光学パラメータで微分することにより、光パルス応答の微分波形を同時に計算する。
【００６２】
本解析法の妥当性を検証するために既に報告されている実験データと解析結果とを比較した結果を図７及び図８に示す。図７に光吸収係数が０．０２（ｍｍ−１）、等価散乱係数が０．９（ｍｍ−１）、厚みがそれぞれ１０（ｍｍ）の散乱体間に、厚みが１０（ｍｍ）の非散乱体を配置した３層不均一散乱体の後方散乱光強度の光源−検出器間隔依存性の解析結果及び実験値を示す。このグラフから本解析結果は、実験値とよく一致することを確認した。図８に後方散乱光パルスの平均遅延時間の光源−検出器間隔依存性の解析結果及び実験値（ファーバンクらによるもの）を示す。なお、表中における「光拡散理論」とは、従来の光拡散方程式を用いた方式による解析値であり、「均一スラブ」とは、非散乱体を含まない均一の散乱体の解析値を示す。
【００６３】
このグラフから本解析結果は、実験値とよく一致することを確認した。以上の結果からも本解析法の妥当性は検証されている。なお、従来の光拡散方程式を用いた方式では、その解が光強度、遅延時間共に光源−検出器間隔が大きくなるに従い実験値との乖離が大きくなり、非散乱体を含む不均一散乱体へは適用ができない。また、図９に非散乱体後方の散乱体吸収係数を変えた場合の後方散乱光パルス計算例を示す。このグラフより、後方散乱光光パルス波形から、非散乱体後方の灰／白質３３に対応する散乱体の光吸収係数変化を容易に同定可能なことが分かる。
【００６４】
続いて、理論値導出部１５が、上記求められた拡散放射流速及び放射発散度等の光伝搬の解析結果を演算することにより、光パルス応答の理論値Ｙ（ｔｎ）を導出する（Ｓ０５、理論値導出ステップ）。
【００６５】
続いて、光吸収係数同定部１６が、上記導出された光パルス応答の理論値Ｙ（ｔｎ）と、光パルス測定装置２０により測定された光パルス応答の実測値Ｉ０（ｔｎ）との自乗残差を以下の式により、波長λ１の場合及び波長λ２の場合それぞれ計算する（Ｓ０６、光吸収係数同定ステップ）。
【数１５】
2005221488.tif 000017
上式でσ２（ｔｎ）は各時間におけるノイズ標準偏差であり、ｎ１及びｎ２は理論値と実測値とをフィッティングする時間の上限及び下限を示す。
【００６６】
上記により計算された自乗残差が例えば１０−４程度以下で、前回計算値との差が無視できる場合、自乗残差は収束している、と見なし、その際の光学パラメータ、即ち光吸収係数μａ及び等価散乱係数μｓ´の値を同定した値とする（Ｓ０７、光吸収係数同定ステップ）。
【００６７】
上記自乗残差が、前回計算値に比べ減少している場合、以下のように光吸収係数同定部１６が光学パラメータの値を再設定し、その値を用いてＳ０４〜Ｓ０７の処理を再度行う。
【００６８】
上記再設定には、式（１５）を各光学パラメータの変化δａｍで偏微分し、各光学パラメータについての停留値条件∂χ２／∂δａｍ＝０を適用し、以下の式により光学パラメータの補正値を求める。
【数１６】
2005221488.tif 000018
式（１６）のＡｋ及びＢｍｋの各要素は次式で与えられる。
【数１７】
2005221488.tif 000019
【数１８】
2005221488.tif 000020
上式で、ｌは未知光学パラメータ数であり、ここでの未知パラメータは、各セルの光吸収係数、等価散乱係数及び光パルス振幅である。但し、光パルス振幅の最大値で規格化した相対振幅を用いる場合は、未知パラメータは、光吸収係数及び等価散乱係数のみとなる。上記の式により求められたδａｍを光学パラメータの補正値として用い、光学パラメータの値を再設定する（Ｓ０８）。
【００６９】
なお、本実施形態では１箇所の測定点についての例を示したが、複数箇所の理論値及び測定値に上記の最小自乗法を適用する場合も同様の処理で光学パラメータを再設定することができる。また、本実施形態では、全てのセルの光学パラメータを未知数として扱っているが、セルを複数統合化したセル群内の光学パラメータを等しいと仮定した簡易な推定も同様の処理で行うことができる。
【００７０】
続いて、成分濃度同定部１７が、上記同定された波長λ１の場合及び波長λ２の場合についての各セルの光吸収係数μａから、以下の式により血中における還元ヘモグロビン濃度ＣＨｂ及び酸化ヘモグロビン濃度ＣＨｂＯ２を求める（Ｓ０９、成分濃度同定ステップ）。
【数１９】
2005221488.tif 000021
【数２０】
2005221488.tif 000022
【数２１】
2005221488.tif 000023
上式でεＨｂλ１，εＨｂλ２，εＨｂＯ２λ１，εＨｂＯ２λ２は、図３に示した、還元ヘモグロビン及び酸化ヘモグロビンの各波長における光吸収係数である。また、μａλ１（ｉ，ｊ，ｋ）及びμａλ２（ｉ，ｊ，ｋ）はＳ０７において、波長λ１，λ２のそれぞれに対して同定した各セルの光吸収係数である。また、成分濃度同定部１７は、上記の値から血流量Ｂ及び血中酸素飽和度Ｓｐを以下の式から求める。
【数２２】
2005221488.tif 000024
【数２３】
2005221488.tif 000025
【００７１】
上記により、本実施形態によれば、散乱体部分を離散化し、光拡散方程式の微分成分を差分成分で置換え、単純な代数計算のみで解析すること等により、解析時間を短くすることができる。また、非散乱体間の光結合を光拡散方程式における等価波源として考慮したこと等により、人頭部３０、即ち非散乱体を含む不均一散乱体に含まれる成分の情報の同定を精度高くすることができる。これにより、光トモグラフィの適用領域が拡大でき、アルツハイマー病の診断と薬餌効果との定量的な評価を始め、人頭部機能解析に広く応用可能となる。もちろん、人頭部以外の不均一散乱体においても広く応用可能となる。
【００７２】
また、本実施形態においては、血中における還元ヘモグロビン及び酸化ヘモグロビンの各成分の濃度を、光吸収係数同定部１６により導出された光吸収係数から求めている。これにより、本実施形態によれば、人頭部３０、即ち非散乱体を含む不均一散乱体に含まれる成分の濃度まで、一連の処理によって求めることができる。なお、求められた光吸収係数を用いた解析を行うことも考えられ、光吸収係数まで求める形態としてもよい。
【００７３】
また、光散乱係数を用いた分析等のために、光吸収係数同定部１６により、光吸収係数と共に光散乱係数を求める形態としてもよい。光散乱係数を求めない場合は、例えば光散乱係数をパラメータではなく固定値として考慮することにより、計算量の削減等を図ることができる。
【００７４】
また、本実施形態においては、離散化を格子状のセルに分割することで行っており、光拡散方程式を時間的に離散化し、ＦＤＴＤ法を用いているが、必ずしもこの方法をとる必要はない。離散化のために設定された各セルにおける光伝搬の解析結果を導出できるように、光拡散方程式の微分成分を差分成分として置換えることができればよい。
【００７５】
また、本実施形態においては、光吸収係数の同定において自乗残差誤差を最小にする基準に基づいて行っているが、必ずしもこの方法をとる必要はない。光吸収係数を同定することができるように、上記理論値と上記実測値の比較をし、参照することができればよい。
【００７６】
光パルス応答を用いた分析には、例えば以下のような手法が考えられ、それらと比較しても上記方法は優れている。
【００７７】
例えば、生体中のフォトン散乱角度、散乱間隔を乱数生成により統計処理し、散乱体の光パルス応答をシミュレーションするモンテカルロ法という手法が考えられる。この方法は、非散乱体を含む不均一散乱体の光パルス応答を計算可能であるが、解析に最低数百時間以上要し、フォトンの到達確率が小さくなるに従い、理論値自体にばらつきが生じる。また、この方法は、生体の光パルス応答の実測値から内部光学パラメータを同定するために必要な生体の光学パラメータによる光パルス応答の微分を計算することが出来ない。
【００７８】
また、散乱体内の光伝搬を、光拡散方程式を用いて解析し、非散乱体を介した光結合をコンピュータグラフィックにおける被写体照度計算公式を用いて評価する、光拡散方程式と被写体照度計算公式とのハイブリッドインテグレーションという方法が考えられる。この方法は、不均一散乱体の光伝搬解析に適用するために、光拡散方程式の有限要素法（Finite Element Method：ＦＥＭ）による定式化を行い、これに被写体照度計算公式を組み込んでいる。しかし、この方法は、等方散乱光を光結合に組み込んでいないため、光伝搬特性推定精度が悪い。また、この方法は周波数領域で数値解析しており、光パルス応答を算出するためには、フーリエ逆変換を行う必要があり、多大な時間を必要とする。
【００７９】
また、光輸送方程式のＦＤＴＤ解析という方法が考えられる。光輸送方程式は散乱体中の光伝搬に拡散方程式を用いていないので、非散乱体を含む不均一散乱体の光伝搬解析では、光拡散方程式より精度はよい。しかし、この方法では、光輸送方程式をＦＤＴＤ解析するために、光輸送関数の未知関数である光強度を角度方向にも離散化する必要があるため、膨大な計算時間を必要とする。しかも、非散乱体中での光輸送方程式の成立は保証されていない。この方法による非散乱体を含む不均一散乱体の光パルス応答の計算結果は報告されておらず、また実験結果等の妥当性の検証はなされていない。
【００８０】
更に、上記ハイブリッドインテグレーション、及び光輸送方程式のＦＤＴＤ解析による生体内部の光学パラメータ推定法は、平均光強度及び平均遅延時間を用いている。この場合、光吸収係数及び等価散乱係数の推定に最低２組の光源及び検出器を用い、これらの位置を変えたデータが必要となる。このため、上記の方法では多数の光源及び検出器を必要とし、かつ精度が悪い。本実施形態においては、光源及び検出器は１組でよい。
【００８１】
本実施形態では、波長は８００ｎｍ程度の光を使用しているが、セルの格子サイズΔｚは０．５ｍｍとしている。波長に対して格子サイズが大きい値となっているが、これは、光拡散方程式は光の位相情報を消去し、光強度のみを取り扱っているためである。
【００８２】
光拡散方程式のＦＤＴＤ解析において、時間間隔ΔｔとΔｚとは次式で与えられる。
【数２４】
2005221488.tif 000026
上式において、ｃは散乱体中の光速、κは定数である。κは小さいほど、一定の格子サイズに対してΔｔを大きく取れるというメリットが出てくる。図１０にκと散乱光パルス波形の関係を示す。通常光トモグラフィは相対振幅が１０−３〜１０−５程度の振幅を使用するため、κ＝１で十分な精度が確保できる。
【００８３】
図１１に、κ＝１の場合についてΔｚ＝０．１２５ｍｍを基準として、Δｚを変えたときのパルス波形を示す。同図から、光パルス立ち上がりを除いてΔｚ＝０．１２５ｍｍからの偏差は、Δｚ＝０．２５ｍｍに対し０．４％以下、Δｚ＝０．５ｍｍに対し１．４％以下、Δｚ＝１．０ｍｍに対して４．２％以下に収まっていることが分かる。光パルス立ち上がり部は、実験データにおいても高速パルス変調した場合は再現性に欠けるので、光トモグラフィでは、１０−３程度の立ち上がり部はデータフィッティングに使用しないのが妥当である。Δｚ＝１．０ｍｍは、散乱体内部の重要な情報を含むパルステール部での数値誤差が４．２％あり、許容できないと考えられる。
【００８４】
上記の結果から、本実施形態でのΔｚは０．５ｍｍが適切である。
【００８５】
上述した一連の処理を成分情報同定装置１０に実行させるための成分情報同定プログラムについて説明する。図１２に示すように、成分情報同定プログラム４１は、成分情報同定装置１０内の記録媒体４０に形成されたプログラム格納領域４０ａ内に格納されている。
【００８６】
成分情報同定プログラム４１は、図１２に示すように、設定モジュール４１ａと、離散化モジュール４１ｂと、解析結果導出モジュール４１ｃと、等価波源組込モジュール４１ｄと、理論値導出モジュール４１ｅと、光吸収係数同定モジュール４１ｆと、成分濃度同定モジュール４１ｇとを備えて構成される。
【００８７】
設定モジュール４１ａを実行させることによって実現する機能は、成分情報同定装置１０の有する設定部１１の機能と同様である。同様に、離散化モジュール４１ｂは離散化部１２の機能に、解析結果導出モジュール４１ｃは解析結果導出部１３の機能に、等価波源組込モジュール４１ｄは等価波源組込部１４の機能に、理論値導出モジュール４１ｅは理論値導出部１５の機能に、光吸収係数同定モジュール４１ｆは光吸収係数同定部１６の機能に、成分濃度同定モジュール４１ｇは成分濃度同定部１７の機能にそれぞれ対応している。
【００８８】
成分情報同定プログラム４１は、その一部若しくは全部が、通信回線等の伝送媒体を介して伝送され、他の機器により受信されて記録（インストールを含む）される構成としてもよい。
【００８９】
また、上記のプログラムは、コンピュータ読取可能な記録媒体に記録して頒布することが可能である。このような記録媒体には、例えばハードディスク及びフレキシブルディスク等の磁気媒体、ＣＤ−ＲＯＭ及びＤＶＤ−ＲＯＭ等の光学媒体、フロプティカルディスク等の磁気光学媒体、あるいは、プログラム命令を実行又は格納するように特別に配置された、例えばＲＡＭ、ＲＯＭ、及び半導体不揮発性メモリなどのハードウェアデバイスなどが含まれる。また、このような記録媒体からプログラム等を読み取る記録媒体読取用のドライブ（例えばフレキシブルディスクドライブ等）をＣＰＵ１０１に対して接続しておいてもよい。
【図面の簡単な説明】
【００９０】
【図１】本発明の実施形態における成分情報を同定するシステムの構成図である。
【図２】酸化ヘモグロビン及び還元ヘモグロビンの光吸収係数を示すグラフである。
【図３】実施形態における成分情報同定装置の機能構成図である。
【図４】実施形態における成分情報同定装置のハードウェア構成を示すブロック図である。
【図５】実施形態における成分情報同定装置での処理を示すフローチャートである。
【図６】非散乱体部分を伝搬する光を示した図である。
【図７】３層不均一散乱体の後方散乱光強度の光源−検出器間隔依存性の解析結果及び実験値を示したグラフである。
【図８】３層不均一散乱体の後方散乱光パルスの平均遅延時間の光源−検出器間隔依存性の解析結果及び実験値を示したグラフである。
【図９】非散乱体後方の散乱体吸収係数を変えた場合の後方散乱光パルス計算例を示したグラフである。
【図１０】時間間隔Δｔと散乱光パルス波形の関係を示したグラフである。
【図１１】格子サイズΔｚと散乱光パルス波形の関係を示したグラフである。
【図１２】本発明に係る成分情報同定プログラムの構成を示す図である。
【符号の説明】
【００９１】
１０…成分情報同定装置、１１…設定部、１２…離散化部、１３…解析結果導出部、１４…等価波源組込部、１５…理論値導出部、１６…光吸収係数同定部、１７…成分濃度同定部、２０…光パルス測定装置、２１…パルス発生器、２２ａ，２２ｂ…レーザ光源、２３…光スイッチ、２４…光検出器、２５…波形観測装置、２６ａ，２６ｂ…光ガイド、３０…人頭部、３１…頭蓋骨、３２…脳髄液、３３…灰／白質、４０…記録媒体、４０ａ…プログラム格納領域、４１…成分情報同定プログラム、４１ａ…設定モジュール、４１ｂ…離散化モジュール、４１ｃ…解析結果導出モジュール、４１ｄ…等価波源組込モジュール、４１ｅ…理論値導出モジュール、４１ｆ…光吸収係数同定モジュール、４１ｇ…成分濃度同定モジュール、１０１…ＣＰＵ、１０２…ＲＯＭ、１０３…ＲＡＭ、１０４…外部記憶装置、１０５…入力装置、１０６…表示装置、１０７…通信インターフェース。
【図１】
2005221488.tif 000027
【図２】
2005221488.tif 000028
【図３】
2005221488.tif 000029
【図４】
2005221488.tif 000030
【図５】
2005221488.tif 000031
【図６】
2005221488.tif 000032
【図７】
2005221488.tif 000033
【図８】
2005221488.tif 000034
【図９】
2005221488.tif 000035
【図１０】
2005221488.tif 000036
【図１１】
2005221488.tif 000037
【図１２】
2005221488.tif 000038