2022-01-27脳波実験チュートリアル
====
[説明用powerpoint(この行をクリックしてダウンロード)](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/2022-01-26ERP_tutorial_powerpoint.pdf)

# Matlab
インストールは東大HPから:https://www.u-tokyo.ac.jp/adm/dics/ja/matlabcwl.html
東京大学のメールアドレスで無料でご利用いただけます
<img src="https://i.imgur.com/cRfE4lK.jpg" width="50%" style="display: block; margin: auto;">
MATLABを起動すると、以下のようなウィンドウが現れます。
左に作業中のディレクトリ(ファイルや関数が保存されているフォルダ)、真ん中にコマンドウィンドウ(電卓)、右にワークスペース(変数やデータフレームを一時的に保存)が表示されます。
<img src="https://i.imgur.com/Fffv3Vh.jpg" width="100%" style="display: block; margin: auto;">
試しにコマンドウィンドウに以下のコマンドを実行してください(コピペするのもOK)
```
1+2
x=[1 2 3]
y=[5 3 12]
graph = plot(x,y)
```
Rと同じ感覚で、x,y変数とグラフが作成されていると思います。
<img src="https://i.imgur.com/a2isN8f.jpg" width="50%" style="display: block; margin: auto;">
MatlabとRとの違いは以下にあります:
1. Matlabは多次元の行列作り、計算に強い (脳波みたいな多次元のデータには計算時間が減ります)
<img src="https://jp.mathworks.com/help/matlab/math/nddemo_02_ja_JP.gif" width="50%" style="display: block; margin: auto;">
2. Matlabでは、多層多様な情報(電極名、時間、電位値、被験者名)を一つの構造に保存するのが可能です(cf.Rは基本dataframe)。見た目にはわかりにくいですが、必要になった時にまたそれぞれ取り出して合わせるというデザインになっています。
<img src="http://3.bp.blogspot.com/-pEBIfvHi-34/VUnGI1bfC4I/AAAAAAAAWFs/d1xDAwsYFA4/s1600/%E6%A7%8B%E9%80%A0%E4%BD%93%E3%82%92%E5%88%86%E8%A7%A3.png" width="50%" style="display: block; margin: auto;">
4. Matlabでは、`{}` `[]` `()` `.` などの演算子に多義性を持っています
参考:[MATLAB の演算子と特殊文字](https://jp.mathworks.com/help/matlab/matlab_prog/matlab-operators-and-special-characters.html?s_tid=srchtitle_%255B%255D%20%25E6%25BC%2594%25E7%25AE%2597%25E5%25AD%2590_4)
<img src="https://i.imgur.com/UFplcce.jpg" width="100%" style="display: block; margin: auto;">
6. Matlabでは、argument名とargumentの内容指定が並列に書くことになっています(初心者にはわかりにくいかもしれません)
```
EEG = pop_eegfiltnew(EEG, 'locutoff',0.1,'hicutoff',30,'plotfreqz',1);
```
```
EEGというデータに対して、フィルターをかける
低い周波数を 0.1Hzに設定し
高い周波数を 30Hzに設定し
結果をグラフで表示する(1)
```
8. Rで簡単に照らし合わせることを、matlabではループで検索して処理をかけるようになっています
```
「チャンネル名がCzであるデータを取り出す」
R → チャンネル名というコラムに「Cz」条件と照らし合わせて、あっている行をすべて取り出す
Matlab → 一行目を取り出して、必要な情報を付け足して、チャンネル名と照らし合わせる → 次の行
```
9. Matlabでは、たくさん省略が許されます(数学者の天才ぶり)
```
x=[1 2 3]
y=[5 3 12]
plot(x,y,'-bo')
% '-'が実線、'b'が青色、'o'が円の記号
```
<img src="https://i.imgur.com/TIuUSuy.jpg" width="50%" style="display: block; margin: auto;">
10. Matlabでは、関数を別ファイルで保存することができる。数多くの関数が関数集にまとめ、そしてGUIインターフェースで操作することが可能です(GUIインターフェースで操作するのと、コマンドで操作する両方可能です。しかし、コマンドで処理されたデータはGUIインターフェースに反映されているかと要確認)

# 関数集(toolbox)) EEGLAB + ERPLAB
EEGLABは脳波データ処理に特化した関数集。(UC SanDiegoのすごい先生方)
ERPLABはEEGLABをベースにした、事象関連電位研究に特化した関数集(加算平均に強い!)(UC S DavisのLuck先生)。
EEGLABは以下のHPからダウンロードできます
https://sccn.ucsd.edu/eeglab/download.php
アプリではないので、「インストール」する必要がないです。
その代わりに、好きなところに保存・解凍すればいいです。
(関数を使うたびにそのフォルダをアクセスする必要があるので、アクセスしやすいところに保存してください)
ERPLABは以下のHPからダウンロードできます
https://erpinfo.org/erplab
アプリではないので、「インストール」する必要がないというところはEEGLABと同じですが
EEGLABをベースにしているので、必ずeeglab2021.1 > plugins に保存しないと連動できません。
動作を確認:
1. C:\Users\soundlab\Documents\MATLAB\eeglab2021.1(例) にアクセスしてください
2. 'eeglab.m'という関数ファイルが存在するかを確認
3. コマンドウィンドウに 'eeglab'を実行して(eeglab.mという関数を実行するという意味)、内蔵したGUIインターフェスが起動され、toolbox欄に 'ERPLAB'というタグが付いているかを確認
# 基本の流れ

# トカゲ実験
複合語:キルギストカゲ
1. LHHH-HLL (正しい条件)
2. *LHHH-HHL (複合語アクセント規則を適用したが核の位置が間違ったA)
3. *LHHH-HHH (複合語アクセント規則を適用したが核の位置が間違ったB)
4. *LHHH-LHH (複合語アクセント規則を適用してない違反)
N1: 1st mora(0-363ms), 2nd mora(363-726ms), 3rd mora(726-1089ms), 4th mora(1089-1454ms)
N2: 1st mora(1454-1828ms), 2nd(1828-2203ms), 3nd(2203-2578ms)
## 1. 二人分の生データをダウンロード
直接クリックではなく、右クリックして「名前を付けてリンク先を保存」で、好きなところに保存してください
[0106-1.eeg](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0106-1.eeg)
[0106-1.vhdr](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0106-1.vhdr)
[0106-1.vmrk](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0106-1.vmrk)
[0107-1.eeg](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.eeg)
[0107-1.vhdr](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.vhdr)
[0107-1.vmrk](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.vmrk)
<img src="https://brainvision.com/wp-content/uploads/2021/06/BV-Logo_303x81.gif" width="50%" style="display: block; margin: auto;">
> The BrainVision format consists of three separate files:
>
> A text header file (.vhdr) containing meta data (タイトル類)
> A text marker file (.vmrk) containing information about events in the data (trigger)
> A binary data file (.eeg) containing the voltage values of the EEG (電位値)
## 2. import data (データの入力)
各メーカで取れたデータのフォーマットはそれぞれです。
EEGLABで処理するには、まずEEGLABで通用する.setファイルに変換する必要があります。
File > Import data > Using EEGLAB functions and plugins > From Brain Vis. Rec..vhdr or .ahdr file
<img src="https://i.imgur.com/7vv6869.jpg" width="100%" style="display: block; margin: auto;">
0106-1.vhdrを選択
「Load a Brain Vision Data Exchange format dataset」を空欄のまま > OK
「Dataset info-- pop_newset()」 > (わかりやすい被験者の名前でいい、例えば) S1 > OK

グラフを確認する:
Plot > Channel data (scroll)
どんなノイズが入っているでしょうか?(瞬き、筋電、電極の不具合、汗)

構造体の中身を確認する
ALLEEG
EEG
<img src="https://i.imgur.com/eVU9VxP.jpg" width="30%" style="display: block; margin: auto;">
## 3. 電極位置の指定
生データの中に、電極について名前しか記録されていません。グラフの作成などで位置座標の指定が必要です。
EEG.chanlocsがチャネル名のみ、位置情報が空欄になっている
Edit > Channel locations > (デフォルトで指定すればいい)
C:\Users\soundlab\Documents\MATLAB\eeglab2021.1\plugins\dipfit\standard_BEM\elec\standard_1005.elc
<img src="https://media.springernature.com/m685/springer-static/image/art%3A10.1038%2Fs41598-019-53233-y/MediaObjects/41598_2019_53233_Fig1_HTML.png" width="50%" style="display: block; margin: auto;">
<img src="https://i.imgur.com/iVaFtL4.jpg" width="100%" style="display: block; margin: auto;">
## 4. サンプリング回数を下げる
今回の実験では500Hzでデータを取りましたが(1秒間に500回データを取る)、脳波関連の成分の周波数は1~30Hzにあるので、少しサンプリング回数を下げても十分に波形を再現できます。それで後の計算は速くなります。500Hzから250Hzに下げます。
Tools > Change sampling rate > 250 Hz
<img src="https://media.geeksforgeeks.org/wp-content/uploads/20200907140648/DownSamplingWaves.PNG" width="50%" style="display: block; margin: auto;">
保存する際の選択肢
<img src="https://i.imgur.com/YQU82iW.png" width="100%" style="display: block; margin: auto;">
## 5. offline filter
脳波関連の成分は大体1~30Hz
直流電源、筋電は50~60Hz以上
汗などは0.1Hz以下
減衰方式で特定の周波数帯域を遮断する
Filter the data > Basic FIR filter(new, default)
Lower edge of the frequency pass band (Hz):0.1
Higher edge of the frequency pass band (Hz):30
<img src="https://i.imgur.com/dSX8hmN.jpg" width="100%" style="display: block; margin: auto;">
## 6. re-reference 再基準化
基準電極は脳電位、心拍と筋電の影響を受けない場所につけたいです。よく基準電極とされているのは鼻尖、耳、耳裏の突起や全脳の平均です。(言語研究はよく両耳の平均を使う)
Tools > Re-reference the data > Re-reference to channel(s) > TP9, TP10 (shift + クリック)を選択 > OK
<img src="https://pressrelease.brainproducts.com/wp-content/uploads/Rereferencing_in_Analyzer_1.png" width="100%" style="display: block; margin: auto;">
## 7. create events
ERPLAB上各トライアルの条件名を数値で処理しないといけないので、数値で名前を付け直します
| 条件 | trigger | event |
| -------- | -------- | -------- |
| 正しい | S 2 | 2 |
| 適用違反A | S 4 | 4 |
| 適用違反B | S 8 | 8 |
| 不適用違反 | S 16 | 16 |
| 自然さ判断 | S 32 | - |
ERPLAB > Create EEG EVENTLIST > 警告メッセージをcontinue > Advanced > Event Codeに2, Event Levelに2 > Update Line
2, 4, 8, 16の順で足していきます。
Save list as (例)"tokage_eventlist.txt" > APPLY

Info to be used sa Event Type > Numeric Codes (数値以外を指定するとエラーが出やすい)
<img src="https://i.imgur.com/mI2fuRq.jpg" width="50%" style="display: block; margin: auto;">
グラフを確認すると、S2が2に変わったはずです。
<img src="https://i.imgur.com/RcyBz94.jpg" width="50%" style="display: block; margin: auto;">
## 8. assign bin
binは「切り出しの基準」を指定するという手続きです
基本event listとは似っていますが、組み合わせ技でもっと複雑な比較を可能にしてくれます。
例:
「自然さ判断で正解である条件Aのトライアルをbin 1に指定する」
「自然さ判断で不正解である条件Aのトライアルをbin 2に指定する」
参考:[ERP Bin Operations](https://github.com/lucklab/erplab/wiki/ERP-Bin-Operations)
今回は単純な対応にします(Binの指定は必ず1から連続数値で)
| 条件 | trigger | event |bin|
| -------- | -------- | -------- |-------- |
| 正しい | S 2 | 2 |1 |
| 適用違反A | S 4 | 4 |2 |
| 適用違反B | S 8 | 8 | 3 |
| 不適用違反 | S 16 | 16 |4 |
| 自然さ判断 | S 32 | - | - |
別途指定内容のファイルが必要で、以下のリンクからダウンロードして下さい(右クリック「名前を付けてリンク先を保存」)
[Bin_20210808](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/Bin_20210808.txt)
中身は以下のようです
```
bin1
a
.{2}
bin2
b
.{4}
bin3
c
.{8}
bin4
d
.{16}
```
ERPLAB > Assign bins (BINLISTER) > Load Bin Descriptor File from > (Bin fileを指定) > Run
構造体のEEG.EVENTLIST.eveninfoを確認すると、biniの欄にbinの情報が入っているはずです。

## 9. epoch (binを切り出す)
ERPLAB > EXtract bin-based epochs > **-200 3200**と指定 > run
baselineの指定(-200msから0msの間)

N1: 1st mora(0-363ms), 2nd mora(363-726ms), 3rd mora(726-1089ms), 4th mora(1089-1454ms)
N2: 1st mora(1454-1828ms), 2nd(1828-2203ms), 3nd(2203-2578ms)
自然さ判断の指示が出るまでに1000msの空白

グラフで確認すると、区切りになっているはずです

## 10. artifact
瞬きなど認知機能を関係ない成分を決まった閾値で排除します。
認知機能と関係する成分は普通±数㎶しかないのに、瞬きは簡単に±100㎶を超えてしまいます。
ERPLAB > Artifact detection in epoched data > simple voltage threshold > Voltage limitsを -100 100 と指定 > Accept

グラフで確認すると、超えた電極とトライアルを色でマークされたことがわかります。

(インピーダンスが高いのと、接触不良などのノイズも見えましたが、それは後半で対策を説明します)
## 11. average
いよいよ肝心な平均を取ります
ERPLAB > Compute averaged ERPs
EEG Dateset(s) Index: 今のset順番
Exclude epochs marked during artifact detection
Exclude epochs that contain either "boundary" or invalied events
> RUN
erpnames: S1

ERP構造体の中身を確認してみましょう
<img src="https://i.imgur.com/woqzwtz.jpg" width="50%" style="display: block; margin: auto;">
## 12. グラフを作成
ERPLAB > Plot ERP > Plot ERP waveform
波形図
bin(条件): 1:4
channels: 17:19 (Fz, Cz, Pz)
negative is up
Row(s): 3
Column(s): 1


## 13. 二人目も同じ手順でS2.erpまで算出
[0107-1.eeg](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.eeg)
[0107-1.vhdr](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.vhdr)
[0107-1.vmrk](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/0107-1.vmrk)
## 14. 合わせてgroup erpを算出
読み実験より、被験者間のバラツキも普通に激しいですので、複数の被験者からの加算平均が必要です

ERPLAB > Average across ERP sets (Grand Average)

---
分析の基本はここまでです。
---
# 実験補助の手引き
http://phiz.c.u-tokyo.ac.jp/~t_chen/download/2022-01-27ERP実験手続き.pdf

---
# 問題解決
powerpointに戻ります
参考用コマンド:
(右クリックして、名前をつけてリンク先を保存)
(今現在使っているコマンドです。まとまってないところはあっちこっちです。進化できたらぜひ共有してください~)
[前処理](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/step1_accent2021_1st02_20211109.m)
[Permutation](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/step2_accent2021_1st02_20211109.m)
[LMEのための出力](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/outPutForEprimeMerge_20211109.m)
[R上でLMEをかける](http://phiz.c.u-tokyo.ac.jp/~t_chen/download/2021-11-10_2192_2456ms_LME.R)
## Thanks to
小林由紀先生、中谷裕教先生、矢野雅貴先生、伊藤愛音先生、李佳霖先生、陳柏亨さん、黃竹佑先生、広瀬友紀先生、伊藤たかね先生
### 参考文献
1. 入戸野 宏(2005).心理学のための事象関連電位ガイドブック. 北大路書店. `初心者に超おすすめ`
1. 開 一夫・金山 範明 (編), 河内山 隆紀・松本 敦・宮腰 誠 (著)(2016).『脳波解析入門 EEGLAB と SPM を使いこなす』. 東京: 東京大学出版会. `とりあえず一通りMATLABで走らせるには`
1. Luck, S. J. (2014). An introduction to the event-related potential technique. MIT press. `がっつり理論の基礎を`
1. Luck, S. J., & Kappenman, E. S. (Eds.). (2011). The Oxford handbook of event-related potential components. Oxford university press. `ERP成分と認知機能の研究を概観する`
1. 音成 龍司. (1997). よくわかる脳波判読. てんかん活動, 96-133. `測定中何がノイズなのか`
1. 金子 秀樹. (2002). “脳波の概要と知識”―生体計測屋より見た脳波の活用入門―. 繊維製品消費科学, 43(9), 554-561. `中学生でもわかる脳波計測の基礎の説明`
1. 浅井里佳, 辻本浩章, 建部渉, & 村治雅文. 脳波の旅への誘い脳波の旅への誘い, 1993. `医学生向けですが、わかりやすい言葉で脳波判読の基礎を説明してくれる`
1. 上坂吉則.(2000). MATLAB プログラミング入門. 牧野書店. `数少ない日本語のMATLAB入門書`
1. Wallisch, P., Lusignan, M. E., Benayoun, M. D., Baker, T. I., Dickey, A. S., & Hatsopoulos, N. G. (2014). MATLAB for neuroscientists: an introduction to scientific computing in MATLAB. Academic Press. `神経研究に使える華麗なる技ばっかり!`