JP2019093008A - 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル - Google Patents

脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル Download PDF

Info

Publication number
JP2019093008A
JP2019093008A JP2017227094A JP2017227094A JP2019093008A JP 2019093008 A JP2019093008 A JP 2019093008A JP 2017227094 A JP2017227094 A JP 2017227094A JP 2017227094 A JP2017227094 A JP 2017227094A JP 2019093008 A JP2019093008 A JP 2019093008A
Authority
JP
Japan
Prior art keywords
activity
measurement data
network
feature
estimation model
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.)
Granted
Application number
JP2017227094A
Other languages
English (en)
Other versions
JP2019093008A5 (ja
JP6884344B2 (ja
Inventor
淳一郎 平山
Junichiro Hirayama
淳一郎 平山
小川　剛史
Takashi Ogawa
剛史 小川
大樹 守谷
daiki Moriya
大樹 守谷
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
ATR Advanced Telecommunications Research Institute International
RIKEN Institute of Physical and Chemical Research
Original Assignee
ATR Advanced Telecommunications Research Institute International
RIKEN Institute of Physical and Chemical Research
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by ATR Advanced Telecommunications Research Institute International, RIKEN Institute of Physical and Chemical Research filed Critical ATR Advanced Telecommunications Research Institute International
Priority to JP2017227094A priority Critical patent/JP6884344B2/ja
Priority to EP18880259.9A priority patent/EP3718473A4/en
Priority to CN201880076799.9A priority patent/CN111417342B/zh
Priority to US16/766,965 priority patent/US20210035665A1/en
Priority to PCT/JP2018/041950 priority patent/WO2019102901A1/ja
Publication of JP2019093008A publication Critical patent/JP2019093008A/ja
Publication of JP2019093008A5 publication Critical patent/JP2019093008A5/ja
Application granted granted Critical
Publication of JP6884344B2 publication Critical patent/JP6884344B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H10/00—ICT specially adapted for the handling or processing of patient-related medical or healthcare data
    • G16H10/60—ICT specially adapted for the handling or processing of patient-related medical or healthcare data for patient-specific data, e.g. for electronic patient records
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/0033—Features or image-related aspects of imaging apparatus classified in A61B5/00, e.g. for MRI, optical tomography or impedance tomography apparatus; arrangements of imaging apparatus in a room
    • A61B5/004—Features or image-related aspects of imaging apparatus classified in A61B5/00, e.g. for MRI, optical tomography or impedance tomography apparatus; arrangements of imaging apparatus in a room adapted for image acquisition of a particular organ or body part
    • A61B5/0042—Features or image-related aspects of imaging apparatus classified in A61B5/00, e.g. for MRI, optical tomography or impedance tomography apparatus; arrangements of imaging apparatus in a room adapted for image acquisition of a particular organ or body part for the brain
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves 
    • A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves  involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
    • A61B5/316—Modalities, i.e. specific diagnostic methods
    • A61B5/369—Electroencephalography [EEG]
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
    • A61B5/316—Modalities, i.e. specific diagnostic methods
    • A61B5/369—Electroencephalography [EEG]
    • A61B5/372—Analysis of electroencephalograms
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235—Details of waveform analysis
    • A61B5/7246—Details of waveform analysis using correlation, e.g. template matching or determination of similarity
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H30/00—ICT specially adapted for the handling or processing of medical images
    • G16H30/20—ICT specially adapted for the handling or processing of medical images for handling medical images, e.g. DICOM, HL7 or PACS
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Medical Informatics (AREA)
  • Public Health (AREA)
  • General Health & Medical Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Physics & Mathematics (AREA)
  • Pathology (AREA)
  • Biomedical Technology (AREA)
  • Radiology & Medical Imaging (AREA)
  • Animal Behavior & Ethology (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Biophysics (AREA)
  • Veterinary Medicine (AREA)
  • Primary Health Care (AREA)
  • Epidemiology (AREA)
  • Psychiatry (AREA)
  • Neurology (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Psychology (AREA)
  • Data Mining & Analysis (AREA)
  • Databases & Information Systems (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Physiology (AREA)
  • Signal Processing (AREA)
  • Artificial Intelligence (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)
  • Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)

Abstract

【課題】脳波の信号を用いて様々な脳内ネットワークの活動を比較的高速に推定できる脳内ネットワークの活動推定システムを提供する。【解決手段】被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得する手段と、脳波の計測データを入力データとする特徴推定モデルを構築するとともに、特徴推定モデルを規定するパラメータを決定する手段と、脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出する手段と、機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出する手段と、モジュール毎の特徴量と、脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定する手段とを含む。【選択図】図１

Description

本発明は、脳波の計測データに基づいて脳内ネットワークの活動を推定する技術に関する。
脳活動を非侵襲で計測する手法として、機能的磁気共鳴画像法（functional Magnetic Resonance Imaging：ｆＭＲＩ）や近赤外分光計測（Near Infrared Spectroscopy：ＮＩＲＳ）といった血流計測による方法と、脳電図または脳波図（Electroencephalogram：以下、「ＥＥＧ」とも略称する。）や脳磁図（Magnetoencephalography：以下、「ＭＥＧ」とも略称する。）を計測する電磁場計測による方法との２種類に大別される。説明の便宜上、本明細書では、ＥＥＧ計測およびＭＥＧ計測により計測される信号変化（時間波形）を「脳波」と総称する。
これらの計測手法のうち、ＥＥＧ計測は、可搬性、携帯性、価格、普及可能性などの点において、他の計測手法に比較して有利である。このようなＥＥＧ計測を用いて、脳内のモデルを決定する汎用的な機械学習手法が提案されている（例えば、非特許文献１）。
また、ｆＭＲＩにより計測される脳活動パターン画像を用いて、様々な脳内ネットワークの存在が確認されている（例えば、非特許文献２）。
ｆＭＲＩ計測により取得される脳活動パターン画像を用いることで、脳内ネットワークの活動を直接的に観測できる一方で、ｆＭＲＩ計測に要する装置は大型であり、限られた場所にしか配置できなく、かつ、被験者が計測装置内に留まっていなければならないといった、普及可能性などの問題がある。
そこで、ｆＭＲＩ計測データとＥＥＧ計測データとを用いて、脳内ネットワークの現在の活動を推定するような手法が提案されている（例えば、特許文献１）。特許文献１は、ｆＭＲＩおよびＥＥＧの同時計測データを用いた、ニューロフィードバックシステムを開示する。
米国特許出願公開第２０１４／０１４８６５７号明細書
Jun-ichiro Hirayama, Aapo Hyvarinen, Motoaki Kawanabe, "SPLICE: Fully Tractable Hierarchical Extension of ICA with Pooling", Proceedings of the 34th International Conference on Machine Learning, PMLR 70:1491-1500, 2017. B.T.T. Yeo, F.M. Krienen, J. Sepulcre, M.R. Sabuncu, D. Lashkari, M. Hollinshead, J.L. Roffman, J.W. Smoller, L. Zollei, J.R. Polimeni, et al., "The organization of the human cerebral cortex estimated by functional connectivity", J. Neurophysiol., 106 (2011), pp. 1125-1165
上述したような先行技術においては、脳内ネットワークの活動を推定するための情報を取得するためには、ＥＥＧ計測結果に対する逆問題解析が必要になる。このような逆問題解析は、高精度な結果を得るために脳や頭蓋の形状情報およびセンサの位置情報が必要になるため汎用性が低くなるという課題があり、また、演算処理量が多く、算出される脳内ネットワークの活動の時間粒度が大きくならざるを得ないという課題がある。また、先行技術においては、複数の脳内ネットワークのうち、特定の一部の脳内ネットワークについて活動が推定できるのみであった。
ＥＥＧ計測などの比較的簡便な方法により計測される脳波の信号を用いて、様々な脳内ネットワークの活動を比較的高速に推定できる手法が要望されている。
本発明のある実施の形態に従う脳内ネットワークの活動推定システムは、被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得する取得手段と、脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータ決定手段とを含む。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定システムは、さらに、脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出する第１の特徴量算出手段と、機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出する第２の特徴量算出手段と、第１の特徴量算出手段により算出されるモジュール毎の特徴量と、第２の特徴量算出手段により算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定する決定手段とを含む。
上記の形態において、脳内ネットワークの活動推定システムは、さらに、被験者から計測される脳波の計測データを入力データとして特徴推定モデルに与えるとともに、少なくとも特定の脳内ネットワークの活動を表現する１または複数のモジュールについて特徴量を算出する第３の特徴量算出手段と、第３の特徴量算出手段により算出される特徴量に基づいて特定の脳内ネットワークの活動度合いを示す値を算出する脳活動算出手段とを含んでいてもよい。
上記の形態において、脳内ネットワークの活動推定システムは、さらに、特定の脳内ネットワークの活動度合いを示す値に基づいて、被験者に対するアクションの内容を決定するフィードバック手段を含んでいてもよい。
上記の形態において、複数のエレメントの各々は、入力データの線形結合として算出されてもよい。
上記の形態において、入力データは、所定期間の時系列データであり、第１の特徴量算出手段は、各モジュールの所定期間に亘る特徴量の時系列データを時間方向に集約して、モジュール毎の特徴量として出力してもよい。
上記の形態において、パラメータ決定手段は、脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間（短時間パワーまたは振幅）に強い相関性が生じるように、特徴推定モデルを規定するパラメータを推定してもよい。
上記の形態において、脳内ネットワークの活動推定システムは、さらに、被験者毎に、特徴推定モデルを規定するパラメータと特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す情報とを関連付けて格納するとともに、要求に応じて、特定の被験者に対応するパラメータおよび特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す情報を出力する記憶手段を有していてもよい。
本発明の別の実施の形態に従う脳内ネットワークの活動推定方法は、被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータステップとを含む。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定方法は、脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、算出されるモジュール毎の特徴量と、算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定するステップとを含む。
本発明のさらに別の実施の形態に従う脳内ネットワークの活動推定プログラムは、コンピュータに、被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータステップとを実行させる。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定プログラムは、コンピュータに、さらに、脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、算出されるモジュール毎の特徴量と、算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定するステップとを実行させる。
本発明のさらに別の実施の形態に従う脳内ネットワークの活動推定システムは、被験者から計測された脳波の計測データを取得する取得手段と、被験者に関連付けられて予め決定されたパラメータを用いて、脳波の計測データを入力データとする特徴推定モデルを構築するモデル構築手段とを含む。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定システムは、被験者に関連付けられて予め決定された、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す関連付け情報に基づいて、取得手段により取得された脳波の計測データを入力データとして特徴推定モデルに与えるとともに、特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出する算出手段を含む。特徴推定モデルを規定するパラメータは、脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定される。関連付け情報は、パラメータの決定に用いられた脳波の計測データを入力データとして特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、パラメータの決定に用いられた脳波の計測データと同時に被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される。
本発明のさらに別の実施の形態に従う脳内ネットワークの活動推定方法は、被験者から計測された脳波の計測データを取得するステップと、被験者に関連付けられて予め決定されたパラメータを用いて、脳波の計測データを入力データとする特徴推定モデルを構築するステップとを含む。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定方法は、被験者に関連付けられて予め決定された、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報に基づいて、取得された脳波の計測データを入力データとして特徴推定モデルに与えるとともに、特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出するステップを含む。特徴推定モデルを規定するパラメータは、脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定される。関連付け情報は、パラメータの決定に用いられた脳波の計測データを入力データとして特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、パラメータの決定に用いられた脳波の計測データと同時に被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される。
本発明のさらに別の実施の形態に従う脳内ネットワークの活動推定プログラムは、コンピュータに、被験者から計測された脳波の計測データを取得するステップと、被験者に関連付けられて予め決定されたパラメータを用いて、脳波の計測データを入力データとする特徴推定モデルを構築するステップとを実行させる。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。脳内ネットワークの活動推定プログラムは、コンピュータに、被験者に関連付けられて予め決定された、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報に基づいて、取得された脳波の計測データを入力データとして特徴推定モデルに与えるとともに、特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出するステップを実行させる。特徴推定モデルを規定するパラメータは、脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定される。関連付け情報は、パラメータの決定に用いられた脳波の計測データを入力データとして特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、パラメータの決定に用いられた脳波の計測データと同時に被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される。
本発明のさらに別の実施の形態に従えば、被験者から計測された脳波の計測データを入力することで、当該被験者における特定の脳内ネットワークの活動度合いを示す値を出力する学習済み脳活動推定モデルが提供される。学習済み脳活動推定モデルは、脳波の計測データを入力データとする特徴推定モデルと、関連付け情報とを含む。特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、複数のエレメントの各々は、脳波の計測データの各次元に対応する値の線形結合として規定されている。関連付け情報は、被験者に関連付けられて予め決定された、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定するパラメータを含む。学習済み脳活動推定モデルを構築するための学習工程において、コンピュータに、被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、脳波の計測データを入力データとして特徴推定モデルに与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、特徴推定モデルを規定するパラメータを決定するステップと、決定されたパラメータに従って構築された特徴推定モデルに脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、算出されるモジュール毎の特徴量と、算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、関連付け情報を決定するステップとを実行させる。
本発明のある実施の形態によれば、ＥＥＧ計測などの比較的簡便な方法により計測される脳波の信号を用いて、様々な脳内ネットワークの活動を比較的高速に推定できる。
本実施の形態に従う脳内ネットワークの活動推定方法の概要を示す模式図である。 本実施の形態に従う脳内ネットワークの活動推定システムの概要を示す模式図である。 本実施の形態に従う活動推定システムを構成する処理装置の装置構成を示す模式図である。 本実施の形態に従う活動推定システムにおけるモジュール推定の原理を説明するための模式図である。 図４に示す仮定下におけるパラメータ決定の手順を説明するための模式図である。 本実施の形態に従う活動推定システムにおいて用いられるＥＥＧ計測データを入力とする特徴推定モデルの一例を示す模式図である。 本実施の形態に従う活動推定システムにおけるＥＥＧ計測データから特徴量を算出する全体処理を説明するための模式図である。 公知の脳内ネットワークの分類例を示す模式図である。 本実施の形態に従う特徴推定モデルの各要素の大きさの一例を示す図である。 本実施の形態に従う活動推定システムにおける相関性評価の処理手順を説明するための模式図である。 デフォルトモードネットワーク（ＤＭＮ）に対して相対的に強い相関性を示す計測結果を示す図である。 デフォルトモードネットワーク（ＤＭＮ）の１つサブネットワークに対して相対的に強い相関性を示す計測結果を示す図である。 制御ネットワーク（ＣＯＮ）に対して相対的に強い相関性を示す計測結果を示す図である。 制御ネットワーク（ＣＯＮ）に対して相対的に強い相関性を示す別の計測結果を示す図である。 背側注意ネットワーク（ＤＡＮ）に対して相対的に強い相関性を示す計測結果を示す図である。 背側注意ネットワーク（ＤＡＮ）に対して相対的に強い相関性を示す別の計測結果を示す図である。 他のモジュールについての計測結果を示す図である。 本実施の形態に従う活動推定システムの処理装置に機能構成例を示す模式図である。 本実施の形態に従う脳内ネットワークの活動推定方法の処理手順を示すフローチャートである。 本実施の形態に従うニューロフィードバックの実現例を説明するための模式図である。 本実施の形態に従う活動推定システムをネットワーク化した実現例を説明するための模式図である。 本実施の形態に従う活動推定システムの処理装置に機能構成例を示す模式図である。 本実施の形態に従う活動推定システムにおける学習工程によって取得される脳活動推定器を示す模式図である。 本実施の形態に従う脳内ネットワークの別の活動推定方法の処理手順を示すフローチャートである。
本発明の実施の形態について、図面を参照しながら詳細に説明する。なお、図中の同一または相当部分については、同一符号を付してその説明は繰り返さない。
［Ａ．概要］
まず、本実施の形態に従う脳内ネットワークの活動推定方法の概要について説明する。図１は、本実施の形態に従う脳内ネットワークの活動推定方法の概要を示す模式図である。
図１を参照して、本実施の形態の脳内ネットワークの活動推定方法においては、まず、同一の被験者に対して、脳波と、機能的磁気共鳴画像法（functional Magnetic Resonance Imaging）（以下、「ｆＭＲＩ」とも略称する。）とを同時に計測する（（１）ＥＥＧ／ｆＭＲＩ同時計測）。ＥＥＧ計測により得られるデータ（以下、「ＥＥＧ計測データ」とも称す。）およびｆＭＲＩ計測により得られるデータ（以下、「ｆＭＲＩ計測データ」とも称す。）は、同一の被験者の同一の脳活動を示すものとなる。
本実施の形態におけるＥＥＧ計測データは、被験者の頭部に配置される各センサから計測される信号群であり、各センサが配置される位置近傍の脳活動を反映することになる。各センサはチャネルとも称され、ＥＥＧ計測データは多チャンネルの脳波データに相当する。各センサは、典型的には、一対の電極で構成される。
続いて、ＥＥＧ計測データからモジュール（ネットワーク）を推定する（（２）モジュール推定）。モジュールに含まれる情報は、被験者の脳内に設定される相関性の強い信号源集団の活動を示す特徴量を含む。すなわち、モジュール化することで推定される特徴量は、脳内ネットワークを表現し得る。後述するように、推定される特徴量は非線形性を有する。
各モジュールは、ＥＥＧ計測データの複数チャネルから影響を受けるものとして推定される。基本的には、各モジュールは、脳内の対応する信号源近傍のセンサから得られたＥＥＧ計測データの影響をより強く受けることになる。
さらに、１または複数のモジュールに関連する１または複数の上位ファクタについても推定する。上位ファクタは、複数のモジュールを集合させた活動を反映するものである。
また、ｆＭＲＩ計測データから脳内ネットワークの活動が推定される（（３）脳内ネットワーク活動推定）。
脳内ネットワークとは、安静時ネットワークとも称され、単一の脳領域に属する信号源による、または、複数の空間的に離れた脳領域に属する信号源が協調することによる、特徴的な脳活動パターンの総称である。脳内ネットワークは、主に、安静時ｆＭＲＩを用いて定義される。代表的なものに、Default Mode Network（ＤＭＮ）やDorsal Attention Networkなどがある。例えば、何らかの機能が要求されると、特定の脳領域に属する信号源のみが活性化するだけではなく、複数の脳領域に属する信号源が協調して活性化することもある。このような各機能に応じて１または複数の脳領域に属する信号源の活動を示すパターンが脳内ネットワークである。なお、単一の脳領域に着目する場合には、単に（機能的）関心領域と称されることも多いが、本明細書における「脳内ネットワーク」は、単一の脳領域に着目する場合、および、複数の空間的に離れた脳領域に着目する場合の両方を包含する概念である。
脳内ネットワークの各々は、１または複数の信号源に関連付けられるモジュールと同様のレイヤに位置することになる。
上述したように、ＥＥＧ計測データおよびｆＭＲＩ計測データは、同一の被験者から同時に取得した計測結果であるので、それぞれの計測データは同一の脳内ネットワークを反映するはずである。このような事前知識を利用して、推定されるモジュールおよび／または上位ファクタと、推定される脳内ネットワークとの関連性を評価する（（４）相関性評価）。そして、上位ファクタ、および、その上位ファクタに関連付けられる１または複数のモジュールのうち、推定される脳内ネットワークと相対的に強い相関性を有するものを、当該脳内ネットワークの活動を表現する候補として決定する（（５）脳内ネットワークの活動を表現するモジュールの特定）。
以上のような処理手順によって、脳内ネットワークの活動を表現するモジュールおよび／または上位ファクタを決定できる。説明の便宜上、以下の説明においては、主として、脳内ネットワークの活動を表現するモジュールを推定する処理について説明するが、モジュールだけではなく、脳内ネットワークの活動を表現する上位ファクタを推定するようにしてもよい。
本実施の形態に従う脳内ネットワークの活動推定方法においては、脳の電気活動により生じる電圧の計測値であるＥＥＧに限らず、脳の電気活動により生じる変動磁界の計測値である脳磁図（Magnetoencephalography：ＭＥＧ）を用いることもできる。
上述したように、本明細書において、「脳波」は、ＥＥＧ計測およびＭＥＧ計測により計測される信号変化（時間波形）を総称する用語である。説明の便宜上、以下の説明においては、主として、ＥＥＧ計測データを用いる例について説明する。
［Ｂ．活動推定システムの概要］
次に、本実施の形態に従う脳内ネットワークの活動推定方法を実現するための活動推定システムの概要について説明する。
図２は、本実施の形態に従う脳内ネットワークの活動推定システム１の概要を示す模式図である。図２を参照して、活動推定システム１は、処理装置１００と、ＥＥＧ計測装置２００と、ｆＭＲＩ計測装置３００とを含む。
処理装置１００は、後述するように、ＥＥＧ計測装置２００により計測されたＥＥＧ計測データと、ｆＭＲＩ計測装置３００により計測されたｆＭＲＩ計測データとを受付けて、特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報を出力する。このように、処理装置１００は、被験者から同時に計測された脳波の計測データ（ＥＥＧ計測データ）および機能的磁気共鳴画像法の計測データ（ｆＭＲＩ計測データ）を取得する。
ＥＥＧ計測装置２００は、被験者Ｓの頭部に配置された複数のセンサ２２０に生じる脳波を示す信号（電気信号）を検出する。ＥＥＧ計測装置２００は、マルチプレクサ２０２と、ノイズフィルタ２０４と、Ａ／Ｄ（Analog to Digital）変換器２０６と、記憶部２０８と、インターフェイス２１０とを含む。
マルチプレクサ２０２は、複数のセンサ２２０にそれぞれ接続されるケーブル２２２のうち一組のケーブルを順次選択して、ノイズフィルタ２０４と電気的に接続する。ノイズフィルタ２０４は、高周波カットフィルタなどのノイズを除去するフィルタであり、選択されたチャネルに対応する一組のケーブル間に生じる脳波を示す信号（電気信号）に含まれるノイズ成分を除去する。
Ａ／Ｄ変換器２０６は、ノイズフィルタ２０４から出力される電気信号（アナログ信号）を所定周期毎にサンプリングして、デジタル信号として出力する。記憶部２０８は、Ａ／Ｄ変換器２０６から出力される時系列データ（デジタル信号）を、選択されているチャネルおよびタイミングを示す情報（例えば、時刻またはカウンタ値）と関連付けて逐次格納する。
インターフェイス２１０は、処理装置１００などからのアクセスに応じて、記憶部２０８に格納されている脳波を示す時系列データを処理装置１００へ出力する。
一方、ｆＭＲＩ計測装置３００は、被験者Ｓの脳活動の情報を取得したい領域（以下、「関心領域」とも称す。）に向けて、共鳴周波数の高周波電磁場を印加することで、特定の原子核（例えば、水素原子核）から共鳴により生じる電磁波を検出することで、脳活動を計測する。
ｆＭＲＩ計測装置３００は、磁場印加機構３１０と、受信コイル３０２と、駆動部３２０と、データ処理部３５０とを含む。
磁場印加機構３１０は、被験者Ｓの関心領域に制御された磁場（静磁場および傾斜磁場）を印加するとともに、ＲＦ（Radio Frequency）パルスを照射する。より具体的には、磁場印加機構３１０は、静磁場発生コイル３１２と、傾斜磁場発生コイル３１４と、ＲＦ照射部３１６と、被験者Ｓをボア内に載置する寝台３１８とを含む。
駆動部３２０は、磁場印加機構３１０に接続され、被験者Ｓに印加される磁場、および、ＲＦパルス波の送受信を制御する。より具体的には、駆動部３２０は、静磁場電源３２２と、傾斜磁場電源３２４と、信号送信部３２６と、信号受信部３２８と、寝台駆動部３３０とを含む。
図２においては、被験者Ｓが載置される円筒形状のボアの中心軸をＺ軸と定義し、Ｚ軸に直交する水平方向および鉛直方向をそれぞれにＸ軸およびＹ軸と定義する。
静磁場発生コイル３１２は、Ｚ軸周りに巻回される螺旋コイルから、ボア内にＺ軸方向の静磁場を発生させる。傾斜磁場発生コイル３１４は、ボア内にＸ軸方向、Ｙ軸方向、Ｚ軸方向の傾斜磁場をそれぞれ発生させる、Ｘコイル、Ｙコイル、Ｚコイル（図示していない）を含む。ＲＦ照射部３１６は、制御シーケンスに従って信号送信部３２６から送信される高周波信号に基づいて、被験者Ｓの関心領域にＲＦパルスを照射する。図２には、ＲＦ照射部３１６が磁場印加機構３１０に内蔵されている構成を例示するが、ＲＦ照射部３１６を寝台３１８側に設けてもよいし、ＲＦ照射部３１６と受信コイル３０２とを一体化してもよい。
受信コイル３０２は、被験者Ｓから放出される電磁波（ＮＭＲ信号）を受信し、アナログ信号を出力する。受信コイル３０２から出力されたアナログ信号は、信号受信部３２８において、増幅およびＡ／Ｄ変換された上で、データ処理部３５０へ出力される。受信コイル３０２は、ＮＭＲ信号を高感度で検出できるように、被験者Ｓに近接して配置されることが好ましい。
データ処理部３５０は、駆動部３２０に対する制御シーケンスを設定するとともに、受信コイル３０２で受信されたＮＭＲ信号から脳活動を示す情報として、脳内の活性度を示す複数の脳活動パターン画像を出力する。
データ処理部３５０は、制御部３５１と、入力部３５２と、表示部３５３と、記憶部３５４と、画像処理部３５６と、データ収集部３５７と、インターフェイス３５８とを含む。データ処理部３５０としては、専用のコンピュータであってもよいし、記憶部３５４などに格納された制御プログラムを実行することで、所定の処理を実現する汎用コンピュータであってもよい。
制御部３５１は、駆動部３２０を駆動させる制御シーケンスを発生させるなどの各機能部の動作を制御する。入力部３５２は、図示しない操作者から各種操作や情報入力を受付ける。表示部３５３は、被験者Ｓの関心領域に関する各種画像および各種情報を画面表示する。記憶部３５４は、ｆＭＲＩ計測に係る処理を実行するための制御プログラム、パラメータ、画像データ（３次元モデル像等）、その他の電子データなどを格納する。画像処理部３５６は、検出されたＮＭＲ信号のデータに基づいて、複数の脳活動パターン画像を生成する。インターフェイス３５８は、駆動部３２０との間で各種の信号を遣り取りする。データ収集部３５７は、関心領域に由来する一群のＮＭＲ信号からなるデータを収集する。
図３は、本実施の形態に従う活動推定システム１を構成する処理装置１００の装置構成を示す模式図である。処理装置１００は、典型的には、汎用的なアーキテクチャに従うコンピュータを採用することができる。図３を参照して、処理装置１００は、主たるコンポーネントとして、プロセッサ１０２と、主記憶部１０４と、コントロールインターフェイス１０６と、ネットワークインターフェイス１０８と、入力部１１０と、表示部１１２と、二次記憶部１２０とを含む。
プロセッサ１０２は、ＣＰＵ（Central Processing Unit）やＧＰＵ（Graphical Processing Unit）といった演算処理回路からなり、二次記憶部１２０に格納されている各種プログラムに含まれるコードを指定される順序に実行することで、後述する各種機能を実現する。主記憶部１０４は、ＤＲＡＭ（Dynamic Random Access Memory）などで構成され、プロセッサ１０２で実行されるプログラムのコードやプログラムの実行に必要な各種ワークデータを保持する。
処理装置１００は、通信機能を有しており、この通信機能は、主として、コントロールインターフェイス１０６およびネットワークインターフェイス１０８によって提供される。
コントロールインターフェイス１０６は、ｆＭＲＩ計測装置３００のデータ処理部３５０との間でデータを遣り取りする。ネットワークインターフェイス１０８は、外部装置（例えば、クラウド上のデータサーバ装置など）との間でデータを遣り取りする。コントロールインターフェイス１０６およびネットワークインターフェイス１０８は、有線ＬＡＮ（Local Area Network）、無線ＬＡＮ、ＵＳＢ（Universal Serial Bus）、Ｂｌｕｅｔｏｏｔｈ（登録商標）などの任意の通信コンポーネントで構成される。
入力部１１０は、典型的には、マウスまたはキーボードなどで構成され、ユーザからの操作を受付ける。表示部１１２は、典型的には、ディスプレイなどで構成され、処理装置１００における処理の実行状態や操作に係る各種情報をユーザへ通知する。
二次記憶部１２０は、典型的には、ハードディスクまたはＳＳＤ（Solid State Drive）などで構成され、プロセッサ１０２にて実行される各種プログラム、処理に必要な各種データ、設定値などを保持する。より具体的には、二次記憶部１２０は、ＥＥＧ計測データ１２１と、ｆＭＲＩ計測データ１２２と、モジュール推定プログラム１２３と、脳内ネットワーク活動推定プログラム１２４と、評価プログラム１２５と、生成モデルパラメータ１７０と、関連付け情報１８０とを格納する。
［Ｃ．ＥＥＧ／ｆＭＲＩ同時計測］
まず、図１に示す「（１）ＥＥＧ／ｆＭＲＩ同時計測」について説明する。図２に示す活動推定システム１を利用して、被験者Ｓは、頭部にセンサを装着した状態で、ｆＭＲＩ計測装置３００のボアに載置されて、ＥＥＧ計測およびｆＭＲＩ計測が並列的に実行される。
処理装置１００は、ＥＥＧ計測装置２００およびｆＭＲＩ計測装置３００からの計測データを共通の時刻を基準として互いに対応付けて格納する。このような共通の時刻に基づく計測データの対応付けによって、時間軸を共通とする、ＥＥＧ計測データおよびｆＭＲＩ計測データを取得できる。これらの２種類の計測データを用いて、両者の相関性を評価する。
［Ｄ．モジュール推定］
次に、図１に示す「（２）モジュール推定」について説明する。以下の説明においては、モジュール推定処理の概略を説明するが、より詳細な処理手順などについては、「Jun-ichiro Hirayama, Aapo Hyvarinen, Motoaki Kawanabe, "SPLICE: Fully Tractable Hierarchical Extension of ICA with Pooling", Proceedings of the 34th International Conference on Machine Learning, PMLR 70:1491-1500, 2017.」（非特許文献１）などを参照されたい。
（ｄ１：原理的な説明）
まず、本実施の形態に従う活動推定システム１におけるモジュール推定の原理的な説明を行なう。図４は、本実施の形態に従う活動推定システム１におけるモジュール推定の原理を説明するための模式図である。図５は、図４に示す仮定下におけるパラメータ決定の手順を説明するための模式図である。
図４を参照して、脳内には、脳活動に応じて活性化する複数の信号源（以下、単に「ソース」とも称す。）が存在すると仮定する。ここで、ＥＥＧ計測データの各チャネルの時間波形は、仮定された複数のソースの時間波形の足し合わせとして求められる。すなわち、仮定された複数のソースは、計測されるＥＥＧ計測データの各チャネルの時間波形を決定する要素として仮定される。
適切な周波数帯を見ると、このようなソースの時間信号の短時間パワー（および、包絡線）は、血中酸素濃度に依存する時間信号（以下、「ＢＯＬＤ信号」とも称す。）の短時間パワー（および、包絡線）とそれぞれ強い相関性を有することが報告されている。ソースの時間信号の短時間パワー（および、包絡線）は、ソースの活動度合いを示すソース活動度を示すことになる。このようなソース活動度は、ソースの時間信号の振幅を用いることもできる。
定義より、同一のネットワーク内ではＢＯＬＤ信号間には強い相関性がある。そのため、ＢＯＬＤ信号と各ソース活動度の時間信号との間の強い相関性を仮定すると、同一のネットワークに属するソース活動度の時間信号間にも強い相関性があると仮定できる。
図５に示すように、複数チャネルのＥＥＧ計測データから各ソースの時間波形を算出する空間フィルタを導入する。導入される空間フィルタは、ＥＥＧ計測データと各ソースの時間波形とを対応付ける写像行列となる。本実施の形態においては、一例として、ＥＥＧ計測データと各ソースの時間波形との間は、線形結合であるとしている。
上述のような仮定の下で、複数チャネルのＥＥＧ計測データの確率分布を表現する生成モデル（図１）を採用し、生成モデルに含まれるパラメータを再帰的に変動させることで、算出されるソース活動度の時間波形間で相関性が高くなる組み合わせを探索する。このような生成モデルの詳細については後述する。
図５に示す空間フィルタは、実際には、生成モデルに含まれるパラメータ行列の逆行列の各行として定義され、各ソース活動度の算出は、これらの空間フィルタを用いて行なう。また、相関性の高いソースの組み合わせを「モジュール」と称し、空間フィルタおよびモジュールを用いてＥＥＧ計測データからソースおよびネットワークの活動度を表現し得る特徴量を算出するモデルを「特徴推定モデル１０」と称す。特徴推定モデル１０は、上位ファクタの活動度を表現する特徴量も算出するようにしてもよい。
より具体的には、活動推定システム１は、脳波の計測データ（ＥＥＧ計測データ）を入力データとして与えた場合に、その確率分布をよく再現するように、生成モデルを規定するパラメータを推定する。生成モデルがソース活動度の時間信号間に強い相関性を仮定するため、このような推定は、算出される各エレメントの時間波形間に強い相関性が生じるようにパラメータを決定することと同様の意味をもつ。
特徴推定モデル１０のパラメータは生成モデルのパラメータと共通であるため、このような生成モデルのパラメータ推定処理によって、特徴推定モデル１０を規定するパラメータを決定あるいは推定できる。このように、活動推定システム１は、脳波の計測データ（ＥＥＧ計測データ）を入力データとする特徴推定モデル１０を構築するとともに、特徴推定モデル１０を規定するパラメータを決定するパラメータ決定機能を有している。
本実施の形態に従うパラメータの推定方法によれば、本来的な教師データであるＢＯＬＤ信号が存在しなくても、複数のソースを仮定して再帰処理を適用することで、ＥＥＧ計測データとソースおよびモジュールとの関係を規定する特徴推定モデル１０を決定できる。
（ｄ２：特徴推定モデル１０）
次に、本実施の形態に従う活動推定システム１において用いられる特徴推定モデル１０について説明する。図６は、本実施の形態に従う活動推定システム１において用いられるＥＥＧ計測データを入力とする特徴推定モデル１０の一例を示す模式図である。図６には、特徴推定モデル１０の最小構成の一例を示す。
図６を参照して、本実施の形態においては、第１レイヤおよび第２レイヤからなる特徴推定モデル１０を想定する。第１レイヤには、ＥＥＧ計測データ１２のチャネル数と同数のソース１４が存在すると仮定する。第１レイヤからは、複数のモジュール１６の特徴量が出力される。
特徴推定モデル１０の入力である前処理済みのＥＥＧ計測データ１２をｎ個の入力ベクトルｘｔ（ｔ＝１，２，・・・，ｎ）と総称する。入力ベクトルｘｔの各々には、実数または虚数である、ベクトルｄの要素ｘｉｔが格納される。各入力ベクトルｘｔの要素は、ソースベクトルｓｔの線形結合であるとする。ソースベクトルｓｔの数は、入力ベクトルｘｔと同数であるとする。このように、特徴推定モデル１０は、脳内の信号源を示す複数のエレメント（ソース１４）と、各々が複数のエレメントの少なくとも一部に結合された複数のモジュール１６とを含む。第２レイヤにおいて、複数のモジュール１６と結合される複数の上位ファクタ１８が存在する。
複数のソース１４と入力ベクトルｘｔとの間は、どのような結合形態であってもよい。但し、複数のソース１４と入力ベクトルｘｔとの間を線形結合とすることで、後述するような特徴推定モデル１０を規定するパラメータのパラメータ推定処理を簡素化できる。すなわち、複数のエレメント（ソース１４）の各々は、入力データの線形結合として算出されてもよい。
本実施の形態に従うモジュール推定においては、独立部分空間分析（Independent Subspace Analysis：ＩＳＡ）に類似した手法を採用する。以下の説明においては、入力ベクトルｘｔおよびソースベクトルｓｔは、独立同分布（independently an identically distributed：ＩＩＤ）であるとする。添え字ｔを省略して、入力ベクトルｘとソースベクトルｓとの間の関係は、係数行列Ａ（すなわち、混合行列）を用いて、以下の（１）式のように表わすことができる。
ｘ＝Ａｓ ・・・（１）
係数行列Ａは、正方行列であり、かつ、可逆行列であるので、係数行列Ａの逆行列である行列Ｗ（≡Ａ−１）は、分離行列に相当する。
以下の説明においては、オリジナルデータベクトルの平均値を差し引くことで、入力ベクトルｘおよびソースベクトルｓの平均はいずれもゼロになっているとする。このようなゼロ平均のベクトルを用いても、一般性を失わない。
第１レイヤのソース１４を互いに重複しないｍ個のグループに分割する。ｊ番目のグループ（１≦ｊ≦ｍ）のｄｊのソースからなるベクトルをソースベクトルｓ［ｊ］と示すものとする。ここで、ソースベクトルｓ［ｊ］は、オリジナルデータ空間のうち、係数行列Ａの対応する列を乗じられたｄｊ次元の部分空間を示す。
本実施の形態に係る特徴推定モデル１０において、ｍ個のソースベクトルｓ［ｊ］は、「パワー（powers）」あるいは「エネルギー（energies）」において、互いに依存関係が存在する。このような「パワー」あるいは「エネルギー」は、Ｌ２ノルム（すなわち、二乗和）である||ｓ［ｊ］||２として定量化できる。このような依存関係は、第２レイヤにある別のソースの線形結合を用いることでモデル化される。このようなモデルでは、点次元の非線形性を生じない。具体的には、Ｌ２ノルムは、以下の（２）式のように表わすことができる（但し、ｊ＝１，２，・・・，ｍ）。
||ｓ［ｊ］||２＝Ｆｊ −１（［Ａ’ｓ’］ｊ） ・・・（２）
ここで、単調写像関数であるＦｊは、逆関数Ｆｊ −１によって（非負の）二乗和を実値にマッピングする。Ａ’およびｓ’は、第２レイヤにおける可逆の結合行列およびソースベクトルである。係数行列Ａ’の逆行列である行列Ｗ’（≡Ａ’−１）は、分離行列に相当する。また、添え字ｊは、ベクトルのｊ番目の要素を意味する。
すなわち、第２レイヤにおいては、モジュール１６に対応するベクトルｘ’＝Ａ’ｓ’となる。
説明の便宜上、すべてのソースベクトルｓ［ｊ］について、以下の（３）式のような規格化ベクトルを導入することで、第１レイヤおよび第２レイヤを有する特徴推定モデル１０を規定する。
ｕ［ｊ］＝ｓ［ｊ］／||ｓ［ｊ］|| ・・・（３）
規格化ベクトルｕ［ｊ］は、すべてのｊについて、他のランダム変数から独立して、単位超球面において一様に分布しているとする。
（ｄ３：パラメータ決定）
本実施の形態においては、ＥＥＧ計測データ１２を表す入力ベクトルｘについての確率密度関数（Probability Density Function：ＰＤＦ）に対する最尤推定を適用することで、特徴推定モデル１０を規定するパラメータを決定する。
本実施の形態に従う特徴推定モデル１０において、入力ベクトルｘについての確率密度関数ｐ（ｘ）は、以下の（４）式のように一般化できる。
ここで、Ｈｋは、任意に選択される非ガウス分布の事前想定に対応する固定の関数を意味する。また、（４）式において、Ｗ［ｊ］は行列Ｗのｄｊ列の要素のみからなる部分行列を意味し、ソースベクトルｓ［ｊ］＝Ｗ［ｊ］ｘで算出できる。
上述の確率密度関数ｐ（ｘ）に最尤推定を適用するにあたっては、尤度の最大化に対応する限りにおいて、任意の損失関数を適用できる。本実施の形態においては、損失関数Ｌとして、確率密度関数ｐ（ｘ）を対数化した、Ｌ＝−ｌｎｐ（ｘ）＋ｃｏｎｓｔ．を採用する。
本実施の形態においては、プーリング処理（すなわち、時間についての平均化処理）を採用する。プーリング処理の採用にあたって、予め定められた期間（サンプリング周期の整数倍）を処理単位として、パラメータの推定および特徴量の算出が行なわれる。以下の説明においては、時系列データを所定時間毎に区切った単位を「エポック」とも称す。
より具体的には、本実施の形態に従う２つのレイヤからなる特徴推定モデル１０についての損失関数Ｌ（ｌ）は、以下の（５）式のように表わすことができる。（５）式において、ｌはエポックを特定する変数である。以下の（５）式は、上述の（４）式に各エポックの時間単位を導入することで算出される。
再帰的な処理により、Ｔ個の入力ベクトルｘ１，ｘ２，・・・，ｘＴをエポック単位で損失関数Ｌ（ｌ）に入力することで、各エポックについての損失（ｌ）を全エポック分について算出する。そして、エポック間で算術平均をとった値（すなわち、負の対数尤度）を最小化するようにパラメータ（係数行列Ａ（≡Ｗ−１）および係数行列Ａ’（≡Ｗ’−１））を決定する。
本実施の形態に従う特徴推定モデル１０は、ｍ個の各モジュール内のソース活動度に強い相関性を仮定するものであるため、最尤推定の結果として、仮定に合うような相関性の強いソースの組み合わせが各モジュール内に獲得されることになる。
このように、本実施の形態に従う特徴推定モデル１０を規定するパラメータを推定するための機械学習は、教師無し学習で実現できる。
（ｄ４：モジュール特徴量の算出処理）
次に、上述のような手順によって決定されるパラメータを用いて、ＥＥＧ計測から、脳内ネットワークの活動を表現し得る非線形性の特徴量を算出する処理の全体について説明する。
図７は、本実施の形態に従う活動推定システム１におけるＥＥＧ計測データから特徴量を算出する全体処理を説明するための模式図である。図７を参照して、まず、ＥＥＧ計測により取得されたＥＥＧ計測データに対して、前処理が実行されることで、ソースベクトルが生成される（（１）前処理）。前処理には、アーチファクト除去、周波数フィルタリング、および、主成分分析（Principal Component Analysis：ＰＣＡ）による低次元化などの標準的な処理が含まれる。
入力されるＥＥＧ計測データ（生データ）は、チャネル数分の時系列データを含む。図７に示すベクトル表現として、チャネル数と同数の行と、各時系列データのサンプル数と同数の列とからなる入力ベクトルを用いてもよい。但し、処理負荷を低減する観点から、計測された生データを低次元化してもよい。後述する評価例においては、６３チャネル分のＥＥＧ計測データから４８次元の入力ベクトルを生成している。
このように、特徴推定モデル１０には、チャネル数と同数またはそれ以下のチャネル数分の、所定期間の時系列データからなる入力データが与えられる。このような標準的な処理により得られた入力データに対して、線形空間フィルタが適用される。空間フィルタは、入力データのそれぞれの要素を重み付けして加算することで、ソースベクトルの各要素を算出する。
このように、前処理としては、低次元化のための主成分分析（Principal Component Analysis：ＰＣＡ）を含む線形空間フィルタなどを用いて実現してもよい。
本実施の形態においては、入力ベクトルからソースベクトルへの写像は線形結合であるとしており、入力ベクトルとソースベクトルとは、行数および列数がいずれも同じとしている。但し、これに限られることなく、入力ベクトルとソースベクトルとの間で行列サイズを異ならせてもよい。
続いて、ソースベクトルはｍ個をモジュール（ネットワーク）に分割する（（２）ｍ個のモジュールに分割）。入力データであるＥＥＧ計測データからｍ個のモジュールを生成する処理において、上述したような処理手順に従って決定されたパラメータ（係数行列Ａ（≡Ｗ−１））が用いられる。
ｍ個のモジュールのそれぞれについて、各時間に対応する１または複数の要素についての二乗和が算出される（（３）二乗和を計算）。算出される二乗和の各々は、各サンプル区間における各モジュールのパワーに相当する。
この結果、モジュール毎の二乗和を各時間の要素とする二乗和ベクトルが算出される。すなわち、二乗和を要素とするｍ行×Ｔ列のベクトルが生成される。
さらに、生成された二乗和のベクトルを時間方向にｌ個のエポックに分割する（（４）時間方向にｌ個のエポックに分割）。続いて、各エポック内において二乗和（パワー）をプーリング処理（すなわち、時間についての平均化処理）して得られる結果（平均値）を、各エポックに対応する時間区間の各モジュールの特徴量（パワーまたは活動度）として出力される（（５）各エポック内におけるパワーの平均値を計算）。すなわち、時間Ｔをｌ等分したエポック毎のパワーの平均値がそれぞれ算出される。この結果、各エポックのパワーの平均値を要素とするｍ行×ｌ列のベクトルが特徴量として生成される。なお、ｍ行×ｌ列のベクトルに含まれる要素を対数化したものを特徴量として出力してもよい。
このように、各モジュールの所定期間に亘る特徴量の時系列データを時間方向に集約して、モジュール毎の特徴量として出力される。
図７に示す一連の処理によって、入力データとして与えられるＥＥＧ計測データに含まれる各エポックに相当する時間毎に特徴量が算出される。このような一連の特徴量の算出処理は、ＥＥＧ計測データの長さ分だけ繰り返される（再帰的に実行される）。以上のような処理手順によって算出される特徴量は、図４に示す各モジュールの特徴量に相当する。以下の手順においては、各モジュールの特徴量が用いられる。
このように、活動推定システム１は、脳波の計測データ（ＥＥＧ計測データ）を入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出する特徴量算出機能を有している。
［Ｅ．脳内ネットワークの活動推定］
次に、図１に示す「（３）脳内ネットワークの活動推定」について説明する。
図８は、公知の脳内ネットワークの分類例を示す模式図である。図８に示す分類例は、「B.T.T. Yeo, F.M. Krienen, J. Sepulcre, M.R. Sabuncu, D. Lashkari, M. Hollinshead, J.L. Roffman, J.W. Smoller, L. Zollei, J.R. Polimeni, et al., "The organization of the human cerebral cortex estimated by functional connectivity", J. Neurophysiol., 106 (2011), pp. 1125-1165」（非特許文献２）に記載された内容に沿ったものである。
図８を参照して、ｆＭＲＩの機能的な関心領域（Region Of Interest：ＲＯＩ）（すなわち、脳内ネットワーク）の分類例としては、以下のような７つのネットワーク推定が知られている。
具体的には、安静時の脳内ネットワーク（resting state network：ＲＳＮ）として、（ａ）制御ネットワーク（Control network：ＣＯＮ）、（ｂ）背側注意ネットワーク（Dorsal attention network：ＤＡＮ）、（ｃ）デフォルトモードネットワーク（Default mode network：ＤＭＮ）、（ｄ）大脳辺縁系（ＬＩＭ）、（ｅ）体性運動ネットワーク（Somatomotor network：ＳＭＮ）、（ｆ）腹側注意ネットワーク（Ventral attention network：ＶＡＮ）、および、（ｇ）視覚ネットワーク（Visual network：ＶＩＳ）の７種類からなる。
なお、（ａ）制御ネットワーク（Control network：ＣＯＮ）は前頭頭頂ネットワーク（Frontal Parietal network）と称されることもあり、（ｆ）腹側注意ネットワーク（Ventral attention network：ＶＡＮ）は顕著性ネットワーク（Saliency Network）と称されることもある。
さらに、図８に示す７つのＲＳＮはいくつかのサブネットワークに分割される。（ａ）制御ネットワーク（Control network：ＣＯＮ）については３つのサブネットワークに分割され、（ｃ）デフォルトモードネットワーク（Default mode network：ＤＭＮ）については４つのサブネットワークに分割され、それ以外の各ネットワークについてはいずれも２つのサブネットワークに分割される。この結果、ｆＭＲＩの機能的な関心領域としては、合計で１７個のサブネットワークが設定される。
ｆＭＲＩ計測データ（すなわち、脳活動パターン画像）の各々について、上述したような１７個のＲＯＩのそれぞれに対応する特徴量を用いることとする。このように、本実施の形態においては、出力される脳活動パターン画像に設定される１７個のＲＯＩ毎に画像特徴量が脳内ネットワークの活動を示す情報として出力される。
このように、活動推定システム１は、機能的磁気共鳴画像法の計測データ（ｆＭＲＩ計測データ）に基づいて脳内ネットワーク毎の画像特徴量を算出する特徴量算出機能を有している。
なお、上述の説明においては、非特許文献２に記載されたネットワークおよびサブネットワークについて特徴量を算出する例を示す。しかしながら、現時点においては、脳内ネットワークの定義は完全には確立しておらず、他のネットワークあるいは他の機能的なＲＯＩの定義に従って、特徴量を算出するようにしてもよい。
［Ｆ．相関性評価］
次に、図１に示す「（４）相関性評価」について説明する。この相関性評価においては、特徴推定モデル１０を構成する各モジュール（ネットワーク）と安静時の脳内ネットワーク（ＲＳＮ）との相関性を評価する。
図９は、本実施の形態に従う特徴推定モデル１０の各要素の大きさの一例を示す図である。図９を参照して、入力データであるＥＥＧ計測データ（入力ベクトルｘ）は、４８次元（４８ベクトル）とし、入力ベクトルｘと同次元数（４８次元）のソースベクトルｓを設定する。ソースベクトルｓに関連付けられるモジュール（ネットワーク）およびモジュールに関連付けられる上位レイヤは、いずれも２０次元（２０ベクトル）とする。
なお、特徴推定モデル１０のモジュールの数は、尤度に基づく一セッション抜き取り交差確認（Leave-one-session-out Cross-Validation）によって、事前に決定してもよい。
図１０は、本実施の形態に従う活動推定システム１における相関性評価の処理手順を説明するための模式図である。図１０を参照して、一例として、各チャネルのＥＥＧ計測データのサンプリング周期は２ｍｓｅｃ（サンプリング周波数：５００Ｈｚ）であるとする。このようなＥＥＧ計測データに対して、アーチファクト除去、周波数フィルタリング、および、主成分分析による低次元化といった標準的な前処理を行なった後、上述したような特徴推定モデル１０のパラメータの決定処理を行なう。
例えば、１０サンプリング（２０ｍｓｅｃの時間窓）毎に１エポックを設定する。上述の（５）式に示すような、時間についてのプーリング処理を導入した特徴推定モデル１０に対する損失Ｌｅｐｏｃｈを用いた最尤推定によってパラメータを決定する。そして、決定されたパラメータに従って、特徴推定モデル１０の各モジュール（ネットワーク）の特徴量であるソースベクトルｓのパワー||ｓ［ｊ］||２を算出し、各エポック内で時間方向に平均化したものを各モジュールの特徴量（パワーまたは活動度）として出力する。
上述の図７に示す各エポックあたりのサンプル数を「１０」にすると、最終的なソースベクトルｓのパワー||ｓ［ｊ］||２（時間方向に平均化後）は、２００ｍｓｅｃ毎（サンプリング周波数：５Ｈｚ）に算出されることになる。
算出されるソースベクトルｓのパワーの時系列データに対して、血液動態反応関数（hemodynamic response function：ＨＲＦ）を適用（すなわち、畳み込み積分）した結果を用いて、ＲＳＮとの対応関係を評価する。このようなソースベクトルｓのパワーは、特徴推定モデル１０を構成するモジュール（ネットワーク）毎に算出される。
ｆＭＲＩ計測データは、各被験者について、目を開けた状態で８セッション（１セッションあたり５分以下）を行なった結果を、結合したものを用いた。このとき、ＲＦパルスの照射周期（time to repeat：ＴＲ）は、２．４５秒とし、セッション間で全体パワーを標準化した。また、ＥＥＧ計測においては、８−１２Ｈｚのフィルタ（α波モード）を適用した。
以上のような手順によって準備される、モジュール（ネットワーク）の各々が示す特徴量の時系列データと、ＲＳＮの各々との間の相関性を評価する。すなわち、活動推定システム１は、算出される特徴推定モデル１０のモジュール（ネットワーク）毎の特徴量と、算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、複数のモジュール（ネットワーク）のうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定する決定機能を有している。
［Ｇ．脳内ネットワークの活動を表現するモジュールの特定］
次に、図１に示す「（５）脳内ネットワークの活動を表現するモジュールの特定」について、評価例について説明する。モジュールの特定処理においては、各ＲＳＮについて、特徴推定モデル１０を構成する複数のモジュール（上述の例では、２０個のモジュール）のうち、いずれのモジュールとの相関性が強いかを評価する。言い換えれば、特徴推定モデル１０を構成する複数のモジュールのうち、いずれのモジュールに注目すれば、各ＲＳＮの値を推定できるかを特定する。
このような各ＲＳＮがいずれのモジュールと関連するのかを特定できれば、ｆＭＲＩに比較して計測が容易なＥＥＧ計測データから、各ＲＳＮの活動度を推定でき、例えば、ニューロフィードバックなどへ応用できる。
以下の説明においては、図９に示す特徴推定モデル１０を構成する２０個のモジュール（ネットワーク）に対して、便宜上、＃１〜＃２０の番号を付与して、それぞれを特定することとする。
（ｇ１：計測例１（ＤＭＮ））
上述したような本計測例においては、２つのモジュールがデフォルトモードネットワーク（ＤＭＮ）に対して相対的に強い相関性を示した。
図１１は、デフォルトモードネットワーク（ＤＭＮ）に対して相対的に強い相関性を示す計測結果を示す図である。図１１（ａ）は、デフォルトモードネットワーク（ＤＭＮ）および腹側注意ネットワーク（ＶＡＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１１（ｂ）は、モジュール＃１５についての相関性評価の一例を示す。
図１１（ｂ）に示す相関性評価の横軸は、脳活動パターン画像に設定される１７個のＲＯＩに対応し、縦軸は相関値に対応する。また、相関が有意であるＲＯＩ（本相関性評価においては、並び替え検定および多重比較補正後のｐ値が０．０５より小さかったもの）については「＊」を付加している。また、相関性評価のグラフ中には、算出された相関値を併せて示すことがある。以下の相関性評価の結果例についても同様である。
図１１に示すように、モジュール＃１５は、ＤＭＮに対して相対的に強い正の相関を示すとともに、ＶＡＮに対して相対的に強い負の相関を示していることが分かる。この結果は、モジュール＃１５が示す特徴量に注目することで、ＤＭＮおよびＶＡＮの活性度を推定できることを示している。
図１２は、デフォルトモードネットワーク（ＤＭＮ）の１つサブネットワークに対して相対的に強い相関性を示す計測結果を示す図である。図１２（ａ）は、デフォルトモードネットワーク（ＤＭＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１２（ｂ）は、モジュール＃１８についての相関性評価の一例を示す。
図１２に示すように、モジュール＃１８は、ＤＭＮのうち１つのサブネットワークに対して相対的に強い正の相関を示していることが分かる。この結果は、モジュール＃１８が示す特徴量は、ＤＭＮの活性度の推定に役立つことを示している。
（ｇ２：計測例２（ＣＯＮ））
本計測例においては、３つのモジュールが制御ネットワーク（ＣＯＮ）に対して相対的に強い相関性を示した。
図１３は、制御ネットワーク（ＣＯＮ）に対して相対的に強い相関性を示す計測結果を示す図である。図１３（ａ）は、制御ネットワーク（ＣＯＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１３（ｂ）は、モジュール＃３についての相関性評価の一例を示す。
図１３に示すように、モジュール＃３は、ＣＯＮのうち１つのサブネットワーク（後部帯状皮質＋楔前部（PCC+Precuneus））に対して相対的に強い正の相関を示していることが分かる。
図１４は、制御ネットワーク（ＣＯＮ）に対して相対的に強い相関性を示す別の計測結果を示す図である。図１４（ａ）は、制御ネットワーク（ＣＯＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１４（ｂ）は、モジュール＃２０についての相関性評価の一例を示す。
図１４に示すように、モジュール＃２０についても、ＣＯＮのうち１つのサブネットワーク（後部帯状皮質＋楔前部（PCC+Precuneus））に対して相対的に強い正の相関を示していることが分かる。
また、上述の図１２に示すモジュール＃１８についても、ＣＯＮのうち１つのサブネットワーク（後部帯状皮質＋楔前部（PCC+Precuneus））に対して相対的に強い正の相関を示していることが分かる。
さらに、モジュール＃３，＃１８，＃２０の各々は、背側注意ネットワーク（ＤＡＮ）および体性運動ネットワーク（ＳＭＮ）に対して相対的に強い正の相関を示していることが分かる。
以上のような結果は、モジュール＃３，＃１８，＃２０が示す特徴量は、ＣＯＮの活性度、ＤＡＮの活性度、および、ＳＭＮの活性度の推定に役立つことを示している。
（ｇ３：計測例３（ＤＡＮ））
本計測例においては、３つのモジュールが背側注意ネットワーク（ＤＡＮ）に対して相対的に強い相関性を示した。
図１５は、背側注意ネットワーク（ＤＡＮ）に対して相対的に強い相関性を示す計測結果を示す図である。図１５（ａ）は、背側注意ネットワーク（ＤＡＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１５（ｂ）は、モジュール＃７についての相関性評価の一例を示す。
図１５に示すように、モジュール＃７は、ＤＡＮのうち１つのサブネットワークに対して相対的に強い正の相関を示していることが分かる。
図１６は、背側注意ネットワーク（ＤＡＮ）に対して相対的に強い相関性を示す別の計測結果を示す図である。図１６（ａ）は、背側注意ネットワーク（ＤＡＮ）が活性化している状態に対応する脳活動パターン画像の一例を示す。図１６（ｂ）は、モジュール＃９についての相関性評価の一例を示す。
図１６に示すように、モジュール＃９は、ＤＡＮのうち１つのサブネットワークに対して相対的に強い正の相関を示していることが分かる。
さらに、モジュール＃７および＃９は、体性運動ネットワーク（ＳＭＮ）に対して相対的に強い相関（正または負）を示していることが分かる。
以上のような結果は、モジュール＃７および＃９が示す特徴量は、ＤＡＮの活性度およびＳＭＮの活性度の推定に役立つことを示している。
（ｇ４：計測例４（ＳＭＮ／ＶＩＳ／ＬＩＭ））
本計測例においては、いくつかの他のモジュールが体性運動ネットワーク（ＳＭＮ）、大脳辺縁系（ＬＩＭ）、および／または、視覚ネットワーク（ＶＩＳ）に対して相対的に強い相関性を示した。図１７は、他のモジュールについての計測結果を示す図である。
図１７（ａ）は、モジュール＃６についての相関性評価の一例を示す、図１７（ａ）に示すように、モジュール＃６は、ＳＭＮおよびＶＩＳに対して相対的に強い相関性を示す。
図１７（ｂ）は、モジュール＃８についての相関性評価の一例を示す。図１７（ｂ）に示すように、モジュール＃８は、ＶＩＳに対して相対的に強い相関性を示す。
図１７（ｃ）は、モジュール＃１７についての相関性評価の一例を示す。図１７（ｃ）に示すように、モジュール＃１７は、ＳＭＮおよびＶＩＳに対して相対的に強い相関性を示す。
図１７（ｄ）は、モジュール＃１１についての相関性評価の一例を示す。図１７（ｄ）に示すように、モジュール＃１１は、ＬＩＭに対して相対的に強い相関性を示す。
（ｇ５：計測例５）
上述したモジュール以外のモジュールについては、いずれの脳内ネットワークに対しても、有意な相関性を示さなかった。
［Ｈ．処理装置における機能構成］
次に、本実施の形態に従う活動推定システム１の処理装置１００に実現される機能構成の一例について説明する。
図１８は、本実施の形態に従う活動推定システム１の処理装置１００に機能構成例を示す模式図である。図１８に示す各機能は、処理装置１００のプロセッサ１０２がモジュール推定プログラム１２３と、脳内ネットワーク活動推定プログラム１２４と、評価プログラム１２５とを実行することで実現される。これらのプログラムが、本実施の形態に従う脳内ネットワークの活動推定プログラムに相当する。説明の便宜上、３つのプログラムで構成される実装例を示しているにすぎず、より少ないプログラムで実装してもよいし、より多くのプログラムで実装してもよい。
また、処理装置１００内の１または複数のプロセッサを用いてこれらのプログラムを実行するようにしてもよいし、複数の処理装置が互いに連係してこれらのプログラムを実行するようにしてもよい。後者の場合には、いわゆるクラウドシステムと称される、ネットワーク上に配置された複数のコンピュータを用いるようにしてもよい。さらに、プロセッサがプログラムを実行することで実現する構成（ソフトウェア実装）に代えて、その全部または一部をＦＰＧＡ（Field-Programmable Gate Array）やＡＳＩＣ（Application Specific Integrated Circuit）といったハードワイヤードな構成を用いて実現してもよい。
本実施の形態に従う各プログラムは、ＯＳ（Operating System）が提供する機能を利用するような形態で実装してもよく、そのような場合であっても、本願発明の技術的範囲に含まれ得る。
図１８を参照して、処理装置１００は、前処理モジュール１５０と、パラメータ決定モジュール１５２と、後処理モジュール１５４と、ＲＯＩ特徴量抽出モジュール１５６と、相関値算出モジュール１５８と、評価モジュール１６０とを含む。
これらのモジュールのうち、前処理モジュール１５０と、パラメータ決定モジュール１５２と、後処理モジュール１５４とについては、モジュール推定プログラム１２３（図３）により実現され、ＲＯＩ特徴量抽出モジュール１５６については、脳内ネットワーク活動推定プログラム１２４により実現され、評価プログラム１２５については、評価プログラム１２５により実現される。
前処理モジュール１５０は、ＥＥＧ計測データ１２１から特徴推定モデル１０に入力される入力ベクトルｘを生成する。より具体的には、前処理モジュール１５０は、ＥＥＧ計測データ１２１に含まれる各チャネルの時系列データに対するプーリング処理やフィルタリング処理を行なう。
パラメータ決定モジュール１５２は、特徴推定モデル１０を構成するパラメータ（係数行列または写像関数など）を決定する。具体的には、パラメータ決定モジュール１５２は、前処理モジュール１５０により生成される入力ベクトルｘを特徴推定モデル１０に入力することで算出されるソースベクトルｓの値（パワー／二乗和）をファクタとする損失関数（上述の（５）式など参照）に基づいて、損失関数の値を最小化するパラメータを推定する。
後処理モジュール１５４は、パラメータ決定モジュール１５２により決定されたパラメータに従って、特徴推定モデル１０から出力されるソースベクトルｓの値に対して、後処理を実行する。具体的には、後処理モジュール１５４は、特徴推定モデル１０から出力されるソースベクトルｓの値（二乗和のベクトル）を所定時間毎のエポックに分割し、各エポックについてプーリング処理（時間についての平均化処理）などを実行する。
ＲＯＩ特徴量抽出モジュール１５６は、ｆＭＲＩ計測装置３００から出力されるｆＭＲＩ計測データ（１または複数の脳活動パターン画像）のうち、脳内ネットワークに関連する部分（ＲＯＩ）の画像特徴量を出力する。すなわち、ＲＯＩ特徴量抽出モジュール１５６からは、ＲＳＮ別に画像特徴量が出力される。
相関値算出モジュール１５８は、後処理モジュール１５４から出力されるモジュール別の特徴量とＲＳＮ別の画像特徴量との間の相関値をそれぞれ算出する。相関値算出モジュール１５８からの出力は、モジュールとＲＳＮとの組み合わせ毎の相関値となる。
評価モジュール１６０は、相関値算出モジュール１５８からのモジュールとＲＳＮとの組み合わせ毎の相関値に基づいて、安静時の脳内ネットワーク（ＲＳＮ）の各々と関連性の強いモジュールを特定する情報を関連付け情報１８０として出力する。関連付け情報１８０は、各被験者に関連付けられて予め決定されており、特徴推定モデル１０を構成する複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する。
一般的には、関連付け情報１８０は、特定の脳内ネットワークの活動に対する、特徴推定モデル１０を構成する複数のモジュールのそれぞれからの影響の度合い（関連性）を示す。そのため、各脳内ネットワークの活動は、特徴推定モデル１０を構成する複数のモジュールに現れる特徴量の重み付け線形結合として表わすこともできる。
このように、処理装置１００は、被験者から同時に計測された脳波の計測データ（ＥＥＧ計測データ）および機能的磁気共鳴画像法の計測データ（ｆＭＲＩ計測データ）を取得する。
また、パラメータ決定モジュール１５２により決定されたパラメータ（生成モデルパラメータ１７０）についても出力される。生成モデルパラメータ１７０を用いることで、特徴推定モデル１０を再現できる。
上述したように、本実施の形態に従う活動推定システム１においては、学習工程として、特徴推定モデル１０を構築するための生成モデルパラメータ１７０を決定する処理（上述の［Ｄ．モジュール推定］の内容）、および、特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報１８０を決定する処理（上述の［Ｆ．相関性評価］および［Ｇ．脳内ネットワークの活動を表現するモジュールの特定］の内容）を含むことになる。
［Ｉ．処理手順］
次に、本実施の形態に従う脳内ネットワークの活動推定方法の処理手順について説明する。図１９は、本実施の形態に従う脳内ネットワークの活動推定方法の処理手順を示すフローチャートである。図１９に示す一部のステップは、処理装置１００においてプログラムが実行されることで実現されてもよい。
図１９を参照して、まず、ＥＥＧおよびｆＭＲＩを同時計測して、ＥＥＧ計測データおよびｆＭＲＩ計測データを取得する（ステップＳ１００）。
続いて、処理装置１００は、取得されたＥＥＧ計測データに対して前処理を行なうことで、特徴推定モデル１０に与えられる入力データを生成する（ステップＳ１０２）。続いて、処理装置１００は、初期パラメータを用いて特徴推定モデル１０を決定し（ステップＳ１０４）、当該決定した特徴推定モデル１０に対して入力データを与えることでソースベクトルの値を算出する（ステップＳ１０６）。さらに、処理装置１００は、算出されたソースベクトルの値を含む損失関数の値が収束条件に適合しているか否かを判断する（ステップＳ１０８）。
損失関数の値が収束条件に適合していなければ（ステップＳ１０８においてＮＯ）、処理装置１００は、特徴推定モデル１０の規定するパラメータの値を更新し（ステップＳ１１０）、ステップＳ１０６以下の処理を繰り返す。
損失関数の値が収束条件に適合していれば（ステップＳ１０８においてＹＥＳ）、処理装置１００は、現時点の特徴推定モデル１０を規定するパラメータを最終的な特徴推定モデル１０のパラメータ（生成モデルパラメータ１７０）として出力する（ステップＳ１１２）。
処理装置１００は、ステップＳ１１２において決定されたパラメータにより規定される特徴推定モデル１０から出力されるソースベクトルの値からモジュール毎の特徴量を算出する（ステップＳ１１４）。
ステップＳ１０２〜Ｓ１１４の処理と並行して、あるいは、ステップＳ１１４の後に、処理装置１００は、ｆＭＲＩ計測データのうち、脳内ネットワークに関連する部分（ＲＯＩ）の画像特徴量を出力する（ステップＳ１１６）。
そして、処理装置１００は、ステップＳ１１４において算出されたモジュール毎の特徴量とステップＳ１１６において算出されたＲＯＩ毎の画像特徴量との間の相関値を算出する（ステップＳ１１８）。最終的に、処理装置１００は、ステップＳ１１８において算出された、モジュール別の特徴量とＲＳＮ別の画像特徴量との間の相関値に基づいて、各ＲＳＮの各々と関連性の強いモジュールを示す情報を最終的な関連付け情報１８０として出力する（ステップＳ１２０）。
［Ｊ．応用例］
次に、本実施の形態に従う脳内ネットワークの活動推定システム１により決定された特徴推定モデルおよび推定されたモジュールの情報（関連付け情報１８０）の応用例について説明する。一例として、ニューロフィードバックへの応用が可能である。
ニューロフィードバックは、被験者の現在の脳活動パターンを推定し、その推定された脳活動パターンを目標の脳活動パターンへ近づけるように、当該被験者に対して、何らかのアクションを与える技術を包含する。
ｆＭＲＩを用いた脳活動パターンの計測では、ｆＭＲＩ計測装置３００の大きさが実現への制約となり得るが、本実施の形態に従う脳内ネットワークの活動を表現するモジュールの推定方法を用いることで、ＥＥＧ計測データを用いて、被験者の脳活動パターンを推定できる。すなわち、本実施の形態に従う推定方法は、ニューロフィードバックを普及させるための重要な基盤技術を提供するものである。
図２０は、本実施の形態に従うニューロフィードバックの実現例を説明するための模式図である。図２０を参照して、上述したような処理手順によって決定された生成モデルパラメータに従って特徴推定モデル１０を構築するとともに、被験者から計測されたＥＥＧ計測データを当該構築した特徴推定モデル１０へ与える。
活動推定システム１は、被験者から計測される脳波の計測データ（ＥＥＧ計測データ）を入力データとして特徴推定モデル１０に与えるとともに、少なくとも特定の脳内ネットワークの活動を表現する１または複数のモジュールについて特徴量を算出する。そして、活動推定システム１は、特徴推定モデル１０に従って算出される特徴量に基づいて、特定の脳内ネットワークの活動度合いと、脳内活動ネットワークのターゲットとなる活動度合いとの類似度を示す値を算出する。
一例として、図２０には、特徴推定モデル１０を用いて算出される各モジュールの特徴量のうち、予め推定されている特定のモジュールの特徴量を用いて、デフォルトモードネットワーク（ＤＭＮ）、制御ネットワーク（ＣＯＮ）、および背側注意ネットワーク（ＤＡＮ）などの活動度が推定される例を示す。なお、ＥＥＧ計測データをリアルタイムに計測する場合には、これらの脳内ネットワークの活動もほぼリアルタイムで出力できる。
これらの脳内ネットワークの活動に基づいて、フィードバック内容（ターゲットとの類似度）を決定し、被験者に対して、決定したフィードバック内容を、表示装置６００などの提示装置を介して提示する。このとき、例えば、被験者には、類似度が高くなるほど半径が大きくなる円形の図形を表示するなどして、現在の類似度の程度を把握できるような提示手法を用いる。一方、被験者には、例えば、提示される類似度が所定値以上になった場合には、何らかの報酬（例えば、金銭報酬でもよいし、他の報酬であってもよい）を提供することを事前に通知しておく。このような一連の処理によって、被験者に対するニューロフィードバックを実現できる。
なお、被験者に提示する情報としては、上述したような情報には限定されず、図形情報に限定されるものでもなく、被験者が類似度を認識できる情報であればよい。
このように、活動推定システム１は、特定の脳内ネットワークの活動度合いを示す値に基づいて、被験者に対するアクションの内容を決定するフィードバック機能を有している。
本実施の形態に従う特徴推定モデル１０を用いた脳内ネットワークの活動推定方法においては、特徴推定モデル１０を規定する生成モデルパラメータ１７０および特徴推定モデル１０において各ＲＳＮの活動度を示すモジュールの情報（関連付け情報１８０）があれば、リアルタイムで脳内ネットワークの活動を推定できる。このような利点を利用して、例えば、専用の設備を用いて、ＥＥＧとｆＭＲＩとの同時計測を一旦行った後、任意の場所で、ニューロフィードバックに基づく治療を受けることもできる。
図２１は、本実施の形態に従う活動推定システム１をネットワーク化した実現例を説明するための模式図である。図２１を参照して、例えば、専用の計測ステーションにおいて、各被験者に対して、ＥＥＧとｆＭＲＩとの同時計測を行なって、処理装置１００が各被験者についての特徴推定モデル１０を構築する生成モデルパラメータ１７０および対応する関連付け情報１８０を出力する。このような生成モデルパラメータ１７０および関連付け情報１８０は、計測ステーションからサーバ装置４００へ送信される。
サーバ装置４００においては、被験者毎の生成モデルパラメータ１７０および関連付け情報１８０を格納する被験者データ４０２を保持する。
計測ステーションの他、１または複数の治療所のうち各被験者が希望する治療所から、サーバ装置４００へアクセスして、各被験者に対応する生成モデルパラメータ１７０および関連付け情報１８０を取得する。そして、各治療所においては、後述するような処理装置５００が配置されており、取得された生成モデルパラメータ１７０および関連付け情報１８０に基づいて、図２０に示すようなニューロフィードバックを用いた治療が実施される。
このように、サーバ装置４００は、被験者毎に、特徴推定モデル１０を規定する生成モデルパラメータ１７０と特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す情報とを関連付けて格納する（被験者データ４０２）。また、サーバ装置４００は、要求に応じて、特定の被験者に対応する生成モデルパラメータ１７０および特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報１８０を出力する。
図２１に示すようなネットワークシステムを用いることで、ニューロフィードバックを用いた治療（疾患治療）の普及をより促進できる。
次に、図２１に示す本実施の形態に従う活動推定システム１の処理装置５００に実現される機能構成の一例について説明する。処理装置５００の装置構成は、上述の図３に示す処理装置１００の装置構成と同様であるので詳細な説明は繰り返さない。
図２２は、本実施の形態に従う活動推定システム１の処理装置５００に機能構成例を示す模式図である。図２２に示す各機能は、処理装置５００のプロセッサが脳内ネットワーク活動推定プログラム（図３に示す脳内ネットワーク活動推定プログラム１２４と同様）を実行することで実現される。
処理装置５００内の１または複数のプロセッサを用いてこれらのプログラムを実行するようにしてもよいし、複数の処理装置が互いに連係してこれらのプログラムを実行するようにしてもよい。後者の場合には、いわゆるクラウドシステムと称される、ネットワーク上に配置された複数のコンピュータを用いるようにしてもよい。さらに、プロセッサがプログラムを実行することで実現する構成（ソフトウェア実装）に代えて、その全部または一部をＦＰＧＡやＡＳＩＣといったハードワイヤードな構成を用いて実現してもよい。
本実施の形態に従う各プログラムは、ＯＳが提供する機能を利用するような形態で実装してもよく、そのような場合であっても、本願発明の技術的範囲に含まれ得る。
図２２を参照して、処理装置５００は、入力モジュール５６０と、前処理モジュール５５０と、モデル構築モジュール５５２と、特徴推定モデル１０と、後処理モジュール５５４と、脳活動算出モジュール５５８とを含む。
入力モジュール５６０は、被験者から計測された脳波の計測データであるＥＥＧ計測データ１２１を取得する。入力モジュール５６０は、ＥＥＧ計測装置と直接的に接続されてもよいし、ＥＥＧ計測装置の一部として実装されてもよい。あるいは、入力モジュール５６０は、ＥＥＧ計測装置から任意の媒体を介して、ＥＥＧ計測データ１２１を受入れるようにしてもよい。
前処理モジュール５５０は、ＥＥＧ計測データ１２１から特徴推定モデル１０に入力される入力ベクトルｘを生成する。基本的な機能は、図１８に示す前処理モジュール１５０と同様であるので、詳細な説明は繰り返さない。
モデル構築モジュール５５２は、事前のパラメータ推定処理によって決定された、各被験者に関連付けられて予め決定された生成モデルパラメータ１７０を用いて、脳波の計測データを入力データとする特徴推定モデル１０を構築する。すなわち、モデル構築モジュール５５２は、予め用意されたネットワーク構造の各要素に対して、対応するパラメータを設定する。
後処理モジュール５５４は、前処理モジュール５５０からの入力ベクトルが特徴推定モデル１０へ入力されることで算出される、各ソースベクトルｓの値に対して、後処理を実行する。具体的には、後処理モジュール５５４は、特徴推定モデル１０から出力されるソースベクトルｓの値（二乗和のベクトル）を所定時間毎のエポックに分割し、各エポックについてプーリング処理（時間についての平均化処理）などを実行する。
脳活動算出モジュール５５８は、関連付け情報１８０に基づいて、後処理モジュール５５４から出力される各モジュールの特徴量から、各脳内ネットワークの活動度合いを示す値をそれぞれ出力する。関連付け情報１８０は、複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定するための情報（１または複数のパラメータ）である。脳活動算出モジュール５５８から出力される、脳内ネットワークの活動度合いを示す値は、上述したようなニューロフィードバックなどに利用されてもよい。
このように、前処理モジュール５５０、後処理モジュール５５４および脳活動算出モジュール５５８は、被験者に関連付けられて予め決定された関連付け情報１８０に基づいて、当該被験者から取得された脳波の計測データ（ＥＥＧ計測データ１２１）を入力データとして特徴推定モデル１０に与えるとともに、特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出する。
図２３は、本実施の形態に従う活動推定システム１における学習工程によって取得される脳活動推定器２０を示す模式図である。図２３を参照して、脳活動推定器２０は、被験者から計測された脳波の計測データが入力されることで、当該被験者における特定の脳内ネットワークの活動度合いを示す値を出力する。
脳活動推定器２０は、被験者から取得された脳波の計測データを用いた機械学習により決定される特徴推定モデル１０（第１のモデル）と、被験者から取得されたｆＭＲＩ計測データを用いた機械学習により決定される関連付け情報１８０（第２のモデル）とから構成される、学習済み脳活動推定モデルである。
脳活動推定器２０に含まれる特徴推定モデル１０は、脳内の信号源を示す複数のエレメント（ソース）と、各々が複数のエレメントの少なくとも一部に結合された複数のモジュールとを含む。複数のエレメントの各々は、脳波の計測データの各次元に対応する値の線形結合として規定されている。
脳活動推定器２０に含まれる関連付け情報１８０は、被験者に関連付けられて予め決定される。関連付け情報１８０は、特徴推定モデル１０を構成する複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する。なお、特徴推定モデル１０が上位ファクタを推定する場合には、関連付け情報１８０として、特定の脳内ネットワークの活動を表現する１または複数の上位ファクタを特定する１または複数のパラメータを用いてもよい。
このような、脳活動推定器２０を採用することで、ＥＥＧ計測などの比較的簡便な方法により計測される脳波の信号を用いて、様々な脳内ネットワークの活動を比較的高速に推定できる。
次に、本実施の形態に従う活動推定システム１の処理装置５００における脳内ネットワークの活動推定方法の処理手順について説明する。図２４は、本実施の形態に従う脳内ネットワークの別の活動推定方法の処理手順を示すフローチャートである。図２４に示す一部のステップは、処理装置５００においてプログラムが実行されることで実現されてもよい。
図２４を参照して、まず、処理装置５００は、任意の被験者についてのＥＥＧ計測データを取得する（ステップＳ２００）。続いて、処理装置５００は、当該被験者について予め決定された生成モデルパラメータ１７０および関連付け情報１８０を取得し（ステップＳ２０２）、特徴推定モデル１０を構築する（ステップＳ２０４）とともに、特徴推定モデル１０の各モジュールが示す特徴量から各脳内ネットワークの活動度合いを示す関係式（モデル）を構築する（ステップＳ２０６）。
続いて、処理装置５００は、取得されたＥＥＧ計測データに対して前処理を行なうことで、特徴推定モデル１０に与えられる所定期間の入力データ（時系列データ）を生成し（ステップＳ２０８）、生成された入力データを特徴推定モデル１０に与えることで現れる各モジュールの特徴量を算出する（ステップＳ２１０）。そして、処理装置５００は、各モジュールが示す特徴量から各脳内ネットワークの活動度合いを示す値をそれぞれ算出する（ステップＳ２１２）。
処理装置５００は、ステップＳ２００において取得したＥＥＧ計測データのうち未処理のものが残っているか否かを判断し（ステップＳ２１４）、未処理のデータが残っている場合（ステップＳ２１４においてＹＥＳ）には、ステップＳ２０８以下の処理を繰り返す。一方、未処理のデータが残っていない場合（ステップＳ２１４においてＮＯ）には、処理装置５００は、処理を終了する。
なお、ステップＳ２１２において順次算出される脳内ネットワークの活動度合いを示す値を用いた任意の処理を付加することができる。
［Ｋ．利点］
本実施の形態に従う脳内ネットワークの活動推定システムによれば、多チャンネル計測された脳波の信号を解析して、脳内の複数の信号源を前提とした脳内ネットワークの活動度を定量化できる。本実施の形態においては、特徴推定モデルにおいて、脳内の信号源の位置を明示的に決定する必要がないので、原理的に、脳や頭蓋の形状計測や逆問題解析などが不要となる。すなわち、汎用性が高く、かつ、演算処理量を減らして高速な処理を実現できる。
本実施の形態に従う脳内ネットワークの活動推定システムを用いた評価実験においては、先行技術において知られている脳内ネットワークに対して、比較的強い相関値を示すものが存在していることが示された。
本実施の形態に従う脳内ネットワークの活動推定システムによれば、脳内ネットワークの活動度合いを定量化して取得できるので、精神疾患や発達障害などの定量指標（バイオマーカ）の開発も容易化する。
今回開示された実施の形態は、すべての点で例示であって制限的なものではないと考えられるべきである。本発明の範囲は、上記した実施の形態の説明ではなくて特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内でのすべての変更が含まれることが意図される。
１ 活動推定システム、１０ 特徴推定モデル、１２ ＥＥＧ計測データ、１４ ソース、１６ モジュール、１８ 上位ファクタ、２０ 脳活動推定器、１００，５００ 処理装置、１０２ プロセッサ、１０４ 主記憶部、１０６ コントロールインターフェイス、１０８ ネットワークインターフェイス、１１０，３５２ 入力部、１１２，３５３ 表示部、１２０ 二次記憶部、１２１ ＥＥＧ計測データ、１２２ ｆＭＲＩ計測データ、１２３ モジュール推定プログラム、１２４ 脳内ネットワーク活動推定プログラム、１２５ 評価プログラム、１５０，５５０ 前処理モジュール、１５２ パラメータ決定モジュール、１５４，５５４ 後処理モジュール、１５６ 特徴量抽出モジュール、１５８ 相関値算出モジュール、１６０ 評価モジュール、１７０ 生成モデルパラメータ、１８０ 関連付け情報、２００ ＥＥＧ計測装置、２０２ マルチプレクサ、２０４ ノイズフィルタ、２０６ Ａ／Ｄ変換器、２０８，３５４ 記憶部、２１０，３５８ インターフェイス、２２０ センサ、２２２ ケーブル、３００ ｆＭＲＩ計測装置、３０２ 受信コイル、３１０ 磁場印加機構、３１２ 静磁場発生コイル、３１４ 傾斜磁場発生コイル、３１６ 照射部、３１８ 寝台、３２０ 駆動部、３２２ 静磁場電源、３２４ 傾斜磁場電源、３２６ 信号送信部、３２８ 信号受信部、３３０ 寝台駆動部、３５０ データ処理部、３５１ 制御部、３５６ 画像処理部、３５７ データ収集部、４００ サーバ装置、４０２ 被験者データ、５５２ モデル構築モジュール、５５８ 脳活動算出モジュール、５６０ 入力モジュール、６００ 表示装置、Ｓ 被験者。

Claims (13)

  1. 被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得する取得手段と、
    前記脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータ決定手段とを備え、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出する第１の特徴量算出手段と、
    前記機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出する第２の特徴量算出手段と、
    前記第１の特徴量算出手段により算出されるモジュール毎の特徴量と、前記第２の特徴量算出手段により算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定する決定手段とを備える、脳内ネットワークの活動推定システム。
  2. 被験者から計測される脳波の計測データを入力データとして前記特徴推定モデルに与えるとともに、少なくとも前記特定の脳内ネットワークの活動を表現する１または複数のモジュールについて特徴量を算出する第３の特徴量算出手段と、
    前記第３の特徴量算出手段により算出される特徴量に基づいて前記特定の脳内ネットワークの活動度合いを示す値を算出する脳活動算出手段とをさらに備える、請求項１に記載の脳内ネットワークの活動推定システム。
  3. 前記特定の脳内ネットワークの活動度合いを示す値に基づいて、前記被験者に対するアクションの内容を決定するフィードバック手段をさらに備える、請求項２に記載の脳内ネットワークの活動推定システム。
  4. 前記複数のエレメントの各々は、前記入力データの線形結合として算出される、請求項１〜３のいずれか１項に記載の脳内ネットワークの活動推定システム。
  5. 前記入力データは、所定期間の時系列データであり、
    前記第１の特徴量算出手段は、各モジュールの前記所定期間に亘る特徴量の時系列データを時間方向に集約して、モジュール毎の特徴量として出力する、請求項１〜４のいずれか１項に記載の脳内ネットワークの活動推定システム。
  6. 前記パラメータ決定手段は、前記脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、前記特徴推定モデルを規定するパラメータを決定する、請求項１〜５のいずれか１項に記載の脳内ネットワークの活動推定システム。
  7. 被験者毎に、前記特徴推定モデルを規定するパラメータと特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す情報とを関連付けて格納するとともに、要求に応じて、特定の被験者に対応するパラメータおよび特定の脳内ネットワークの活動を表現する１または複数のモジュールを示す情報を出力する記憶手段をさらに備える、請求項１〜６のいずれか１項に記載の脳内ネットワークの活動推定システム。
  8. 被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、
    前記脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータステップとを備え、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、
    前記機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、
    前記算出されるモジュール毎の特徴量と、前記算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定するステップとを備える、脳内ネットワークの活動推定方法。
  9. 脳内ネットワークの活動推定プログラムであって、コンピュータに、
    被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、
    前記脳波の計測データを入力データとする特徴推定モデルを構築するとともに、当該特徴推定モデルを規定するパラメータを決定するパラメータステップとを実行させ、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、
    前記機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、
    前記算出されるモジュール毎の特徴量と、前記算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを決定するステップとを実行させる、脳内ネットワークの活動推定プログラム。
  10. 被験者から計測された脳波の計測データを取得する取得手段と、
    前記被験者に関連付けられて予め決定されたパラメータを用いて、前記脳波の計測データを入力データとする特徴推定モデルを構築するモデル構築手段とを備え、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記被験者に関連付けられて予め決定された、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報に基づいて、前記取得手段により取得された脳波の計測データを入力データとして前記特徴推定モデルに与えるとともに、前記特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出する算出手段を備え、
    前記特徴推定モデルを規定するパラメータは、前記脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定され、
    前記関連付け情報は、前記パラメータの決定に用いられた脳波の計測データを入力データとして前記特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、前記パラメータの決定に用いられた脳波の計測データと同時に前記被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される、脳内ネットワークの活動推定システム。
  11. 被験者から計測された脳波の計測データを取得するステップと、
    前記被験者に関連付けられて予め決定されたパラメータを用いて、前記脳波の計測データを入力データとする特徴推定モデルを構築するステップとを備え、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記被験者に関連付けられて予め決定された、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報に基づいて、前記取得された脳波の計測データを入力データとして前記特徴推定モデルに与えるとともに、前記特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出するステップを備え、
    前記特徴推定モデルを規定するパラメータは、前記脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定され、
    前記関連付け情報は、前記パラメータの決定に用いられた脳波の計測データを入力データとして前記特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、前記パラメータの決定に用いられた脳波の計測データと同時に前記被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される、脳内ネットワークの活動推定方法。
  12. 脳内ネットワークの活動推定プログラムであって、コンピュータに、
    被験者から計測された脳波の計測データを取得するステップと、
    前記被験者に関連付けられて予め決定されたパラメータを用いて、前記脳波の計測データを入力データとする特徴推定モデルを構築するステップとを実行させ、前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、
    前記被験者に関連付けられて予め決定された、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定する関連付け情報に基づいて、前記取得された脳波の計測データを入力データとして前記特徴推定モデルに与えるとともに、前記特定の脳内ネットワークの活動を表現する１または複数のモジュールについての特徴量から当該特定の脳内ネットワークの活動度合いを示す値を算出するステップを実行させ、
    前記特徴推定モデルを規定するパラメータは、前記脳波の計測データを入力データとして与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、決定され、
    前記関連付け情報は、前記パラメータの決定に用いられた脳波の計測データを入力データとして前記特徴推定モデルに与えた場合に生じる各モジュールの出力値に基づいて算出されたモジュール毎の特徴量と、前記パラメータの決定に用いられた脳波の計測データと同時に前記被験者から計測された機能的磁気共鳴画像法の計測データに基づいて算出された脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、決定される、脳内ネットワークの活動推定プログラム。
  13. 被験者から計測された脳波の計測データを入力することで、当該被験者における特定の脳内ネットワークの活動度合いを示す値を出力する学習済み脳活動推定モデルであって、
    前記脳波の計測データを入力データとする特徴推定モデルと、
    関連付け情報とを備え、
    前記特徴推定モデルは、脳内の信号源を示す複数のエレメントと、各々が前記複数のエレメントの少なくとも一部に結合された複数のモジュールとを含み、前記複数のエレメントの各々は、前記脳波の計測データの各次元に対応する値の線形結合として規定されており、
    前記関連付け情報は、前記被験者に関連付けられて予め決定された、前記複数のモジュールのうち特定の脳内ネットワークの活動を表現する１または複数のモジュールを特定するパラメータを含み、
    前記学習済み脳活動推定モデルを構築するための学習工程において、コンピュータに
    被験者から同時に計測された脳波の計測データおよび機能的磁気共鳴画像法の計測データを取得するステップと、
    前記脳波の計測データを入力データとして前記特徴推定モデルに与えた場合に算出される各エレメントの時間波形間に強い相関性が生じるように、前記特徴推定モデルを規定するパラメータを決定するステップと、
    決定されたパラメータに従って構築された特徴推定モデルに前記脳波の計測データを入力データとして与えた場合に算出される各モジュールの出力値に基づいて、モジュール毎の特徴量を算出するステップと、
    前記機能的磁気共鳴画像法の計測データに基づいて脳内ネットワーク毎の画像特徴量を算出するステップと、
    前記算出されるモジュール毎の特徴量と、前記算出される脳内ネットワーク毎の画像特徴量との間の相関性を評価することで、前記関連付け情報を決定するステップとを実行させる、学習済み脳活動推定モデル。
JP2017227094A 2017-11-27 2017-11-27 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル Active JP6884344B2 (ja)

Priority Applications (5)

Application Number Priority Date Filing Date Title
JP2017227094A JP6884344B2 (ja) 2017-11-27 2017-11-27 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル
EP18880259.9A EP3718473A4 (en) 2017-11-27 2018-11-13 INTRACEREBRAL NETWORK ACTIVITY ESTIMATION SYSTEM, INTRACEREBRAL NETWORK ACTIVITY ESTIMATION PROCESS, INTRACEREBRAL NETWORK ACTIVITY ESTIMATION PROGRAM AND CEREBRAL ACTIVITY ESTIMATION MODEL LEARNED
CN201880076799.9A CN111417342B (zh) 2017-11-27 2018-11-13 脑网络的活动估计***、方法和存储介质
US16/766,965 US20210035665A1 (en) 2017-11-27 2018-11-13 Brain network activity estimation system, method of estimating activities of brain network, brain network activity estimation program, and trained brain activity estimation model
PCT/JP2018/041950 WO2019102901A1 (ja) 2017-11-27 2018-11-13 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2017227094A JP6884344B2 (ja) 2017-11-27 2017-11-27 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル

Publications (3)

Publication Number Publication Date
JP2019093008A true JP2019093008A (ja) 2019-06-20
JP2019093008A5 JP2019093008A5 (ja) 2021-01-14
JP6884344B2 JP6884344B2 (ja) 2021-06-09

Family

ID=66631864

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2017227094A Active JP6884344B2 (ja) 2017-11-27 2017-11-27 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル

Country Status (5)

Country Link
US (1) US20210035665A1 (ja)
EP (1) EP3718473A4 (ja)
JP (1) JP6884344B2 (ja)
CN (1) CN111417342B (ja)
WO (1) WO2019102901A1 (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2021133142A (ja) * 2020-02-28 2021-09-13 株式会社日立製作所 医用撮像装置、医用画像処理装置、及び、画像処理プログラム
WO2022004841A1 (ja) 2020-07-02 2022-01-06 株式会社国際電気通信基礎技術研究所 推定システム、推定方法、プログラム、推定モデル、脳活動トレーニング装置、脳活動トレーニング方法、および、脳活動トレーニングプログラム

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11553871B2 (en) * 2019-06-04 2023-01-17 Lab NINE, Inc. System and apparatus for non-invasive measurement of transcranial electrical signals, and method of calibrating and/or using same for various applications
CN112274144A (zh) * 2019-07-22 2021-01-29 苏州布芮恩智能科技有限公司 近红外脑功能成像数据的处理方法、装置和存储介质
CN113017651B (zh) * 2021-03-16 2022-06-21 哈尔滨工业大学 一种情感eeg的脑功能网络分析方法
CN113040789B (zh) * 2021-03-17 2022-09-09 华南理工大学 一种同步EEG-fMRI数据采集中核磁伪迹的在线去除方法
CN115024684A (zh) * 2022-05-23 2022-09-09 天津大学 刺激范式生成***、脑-机接口***、检测方法及装置
CN115359297B (zh) * 2022-08-24 2024-01-26 南京航空航天大学 一种基于高阶脑网络的分类方法、***、电子设备及介质
CN116823813B (zh) * 2023-08-28 2023-12-01 北京航空航天大学 一种基于磁共振影像的脑韧性评估方法、***及电子设备
CN118094277B (zh) * 2024-04-25 2024-07-02 之江实验室 基于动态功能脑网络的大脑视听融合机制探索方法和装置

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005000528A (ja) * 2003-06-13 2005-01-06 National Institute Of Information & Communication Technology 脳機能測定装置及び脳機能測定プログラム
CN103345749A (zh) * 2013-06-27 2013-10-09 中国科学院自动化研究所 一种基于模态融合的大脑网络功能连接偏侧性检测方法
US20140148657A1 (en) * 2011-02-03 2014-05-29 Ramoot At Tel-Aviv University Ltd. Method and system for use in monitoring neural activity in a subject's brain
JP2017192425A (ja) * 2016-04-18 2017-10-26 株式会社アラヤ 認知状態推定システム及び認知状態推定方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9694178B2 (en) * 2013-10-21 2017-07-04 Neuroelectronics Barcelona S.L. Method and a system for optimizing the configuration of multisite transcranial current stimulation and a computer-readable medium
CN102743166B (zh) * 2012-08-09 2013-11-06 西南大学 一种事件相关电位源定位方法
CN107242873B (zh) * 2017-07-05 2020-06-09 成都信息工程大学 一种基于功能磁共振心理生理交互的脑网络构建方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005000528A (ja) * 2003-06-13 2005-01-06 National Institute Of Information & Communication Technology 脳機能測定装置及び脳機能測定プログラム
US20140148657A1 (en) * 2011-02-03 2014-05-29 Ramoot At Tel-Aviv University Ltd. Method and system for use in monitoring neural activity in a subject's brain
CN103345749A (zh) * 2013-06-27 2013-10-09 中国科学院自动化研究所 一种基于模态融合的大脑网络功能连接偏侧性检测方法
JP2017192425A (ja) * 2016-04-18 2017-10-26 株式会社アラヤ 認知状態推定システム及び認知状態推定方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
KARAHAN, ESIN ET AL.: "Tensor Analysis and Fusion of Multimodal Brain Images", PROCEEDINGS OF THE IEEE, vol. 103, no. 9, JPN6019005020, September 2015 (2015-09-01), pages 1531 - 1559, XP011667005, ISSN: 0004478579, DOI: 10.1109/JPROC.2015.2455028 *
MURTA, TERESA ET AL.: "Electrophysiological correlates of the BOLD signal for EEG-informed fMRI", HUMAN BRAIN MAPPING, vol. 36, no. 1, JPN6019005024, 2014, pages 391 - 414, XP055618887, ISSN: 0004478580, DOI: 10.1002/hbm.22623 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2021133142A (ja) * 2020-02-28 2021-09-13 株式会社日立製作所 医用撮像装置、医用画像処理装置、及び、画像処理プログラム
JP7256765B2 (ja) 2020-02-28 2023-04-12 株式会社日立製作所 医用撮像装置、医用画像処理装置、及び、画像処理プログラム
WO2022004841A1 (ja) 2020-07-02 2022-01-06 株式会社国際電気通信基礎技術研究所 推定システム、推定方法、プログラム、推定モデル、脳活動トレーニング装置、脳活動トレーニング方法、および、脳活動トレーニングプログラム

Also Published As

Publication number Publication date
WO2019102901A1 (ja) 2019-05-31
CN111417342B (zh) 2023-06-13
EP3718473A1 (en) 2020-10-07
CN111417342A (zh) 2020-07-14
JP6884344B2 (ja) 2021-06-09
US20210035665A1 (en) 2021-02-04
EP3718473A4 (en) 2021-12-15

Similar Documents

Publication Publication Date Title
WO2019102901A1 (ja) 脳内ネットワークの活動推定システム、脳内ネットワークの活動推定方法、脳内ネットワークの活動推定プログラム、および、学習済み脳活動推定モデル
US20180276822A1 (en) Ultra-dense electrode-based brain imaging system
Auranen et al. Bayesian analysis of the neuromagnetic inverse problem with ℓp-norm priors
Dimitriadis et al. Reliability of static and dynamic network metrics in the resting-state: A MEG-beamformed connectivity analysis
Bridwell et al. Spatiospectral decomposition of multi-subject EEG: evaluating blind source separation algorithms on real and realistic simulated data
Sun et al. Deep neural networks constrained by neural mass models improve electrophysiological source imaging of spatiotemporal brain dynamics
CN110101384B (zh) 用于复杂网络的功能性网络分析***及分析方法
US20160051161A1 (en) Method for locating a brain activity associated with a task
CN114052668B (zh) 一种基于脑磁图数据的脑功能分析方法
JP2018187044A (ja) 感情推定装置、感情推定方法およびコンピュータプログラム
Hansen et al. Unmixing oscillatory brain activity by EEG source localization and empirical mode decomposition
Conte et al. The influence of the head model conductor on the source localization of auditory evoked potentials
WO2022004841A1 (ja) 推定システム、推定方法、プログラム、推定モデル、脳活動トレーニング装置、脳活動トレーニング方法、および、脳活動トレーニングプログラム
Frank et al. Dynamic multiscale modes of resting state brain activity detected by entropy field decomposition
Kaur et al. Recent developments in spatio-temporal EEG source reconstruction techniques
Chang et al. Assessing recurrent interactions in cortical networks: Modeling EEG response to transcranial magnetic stimulation
WO2019028551A1 (en) SYSTEMS AND METHODS FOR QUALITY EVALUATION OF MAGNETIC RESONANCE FUNCTIONAL IMAGING DATA USING FREQUENCY DOMAIN ENTROPY REGULARIZATION
Hu et al. Tensor clustering on outer-product of coefficient and component matrices of independent component analysis for reliable functional magnetic resonance imaging data decomposition
Jiricek et al. Electrical source imaging in freely moving rats: Evaluation of a 12-electrode cortical electroencephalography system
JP2001514031A (ja) 脳および心臓の一次電流の断層撮影システムおよび方法
Vanhatalo et al. Advanced EEG and MRI measurements to study the functional development of the newborn brain
RU2743608C1 (ru) Способ локализации отделов головного мозга
Chintaluri et al. kCSD-python, reliable current source density estimation with quality control
JP6821171B2 (ja) 脳活動状態定量化方法および脳活動状態計測装置
Rashid et al. Analyzing functional magnetic resonance brain images with opencv2

Legal Events

Date Code Title Description
A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20171212

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20171212

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20201120

A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20201120

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20210406

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20210427

R150 Certificate of patent or registration of utility model

Ref document number: 6884344

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250