Browse Source

업데이트 'README.md'

README.md v1.0.0 released
Hio-Been Han 4 years ago
parent
commit
10ae72bb01
1 changed files with 59 additions and 57 deletions
  1. 59 57
      README.md

+ 59 - 57
README.md

@@ -35,13 +35,13 @@ Raw EEG data are saved in EEGLAB dataset format (*.set). Below are the list of f
     
 **d) Example python/matlab scripts**
 
-    [data_description.ipynb (written and tested on Google COLAB - Python 3 environment)
+    [data_description.ipynb, data_description.py (written and tested on Python 3 environment)
     
     
 <br>
 <br>
 # 3. How to get started (Python 3 without _gin_)
-As the data are saved in EEGLAB format, you need to install appropriate module to access the data in Python3 environment. The fastest way would be to use *read_epochs_eeglab()* function in *MNE-python* module. You can download the toolbox from the link below (or use *pip*).
+As the data are saved in EEGLAB format, you need to install appropriate module to access the data in Python3 environment. The fastest way would be to use <code>read_epochs_eeglab()</code> function in *MNE-python* module. You can download the toolbox from the link below (or use <code>pip install mne</code> in terminal shell).
 
 *[MNE-python]* https://martinos.org/mne/stable/index.html
 
@@ -51,11 +51,13 @@ As the data are saved in EEGLAB format, you need to install appropriate module t
 
 ### 1-1. Download dataset and MNE-python module
 
-The dataset has been uploaded on G-Node and can be accessed by git command, by typing *git clone https://gin.g-node.org/hiobeen/Mouse_hdEEG_ASSR_Hwang_et_al*. However, it's currently not functioning because of the large size of each dataset (>100 MB). Instead, you can use _gin_ command or custom function written below to copy dataset into your work environment.
+The dataset has been uploaded on G-Node and can be accessed by git command, by typing <code>git clone https://gin.g-node.org/hiobeen/Mouse_hdEEG_ASSR_Hwang_et_al</code>. However, it's currently not functioning because of the large size of each dataset (>100 MB). Instead, you can use <code>gin</code> command or custom function written below to copy dataset into your work environment. Also, you need to install *MNE-Python* module using *pip* command to load EEGLAB-formatted EEG data.
 
 > Warning: Direct clonning using *git clone git@gin.g-node.org:/hiobeen/Mouse_hdEEG_ASSR_Hwang_et_al.git* may not work because of the large size of each dataset (>100 MB).
 
-Implementation of scripts below will resulted in downloads of dataset files. In addition, you need to install *MNE-Python* module using *pip* command to load EEGLAB-formatted EEG data.
+To download dataset and install MNE-python module into your environment (local machine/COLAB), try running scripts below.
+
+> Note: Through this step-by-step demonstration, we will use data from one animal (#Animal 2). Unnecessary data files will not be downloaded to prevent long download time. To download whole dataset, change this part; <code>dataset_to_download = [2]</code> into <code>dataset_to_download = [1,2,3,4,5,6]</code>.
 
 
 ```python
@@ -115,7 +117,7 @@ system('pip install mne');
 
 ### 1-2. Accessing meta-data table
 
-File *meta.csv* contains the demographic information of 6 mice. Using *read_csv()* of *pandas* module, meta-datat able can be visualized as follow. 
+File *meta.csv* contains the demographic information of 6 mice. Using <code>read_csv()</code> of *pandas* module, meta-datat able can be visualized as follow. 
 
 
 ```python
@@ -128,12 +130,12 @@ meta
 
 ### 1-3. Data loading and dimensionality check
 
-Each _*.fdt_ file is consisted of different number of trials. To load dataset, a function *get_eeg_data()* is defined below. To maintain original dimensionality order (cf. channel-time-trial in EEGLAB of Matlab), *np.moveaxis()* was applied. 
+Each _*.fdt_ file is consisted of different number of trials. To load dataset, a function <code>get_eeg_data()</code> is defined below. To maintain original dimensionality order (cf. channel-time-trial in EEGLAB of Matlab), <code>np.moveaxis()</code> was applied. 
 
 
 ```python
 # Demo 1-3. Data loading and dimensionality check
-from -mne.io import read_epochs_eeglab as loadeeg
+from mne.io import read_epochs_eeglab as loadeeg
 import numpy as np
 def get_eeg_data(dataset_idx, CAL=1e-6):
   f_name = '%s%srawdata/%s'%(dir_origin,dir_dataset,meta.file_name[dataset_idx])
@@ -214,7 +216,8 @@ plt.show();
 ## Part 2. Plotting Event-Related Potentials
 <br>
 ### 2-1. Accessing event info
-Event information is saved in struct-type variable EEG.event and you can access it by typing {EEG.event.type}. Also, their time trace are avilable in 37-th and 38-th channel of EEG.data. As demonstration purpose, light and sound stimuli of 7 different types of event can be extracted and drawn as follow.
+
+Event information is saved in struct-type variable EEG.event and you can access it by typing <code>EEG.event</code>. Also, their time trace are avilable in 37-th and 38-th channel of <code>EEG.data</code>. For demonstration purpose, light and sound stimuli of 7 different types of event can be extracted and drawn as follow.
 
 ```python
 # Demo 2-1. Event profile (sound/light stimuli)
@@ -249,7 +252,7 @@ plt.subplots_adjust(wspace=.3, hspace=.8)
 
 If data is successfully loaded, now you're ready! For data visualization, an example function is provided below. 
 
-The function *plot_multichan()* draws multi-channel time series data, by taking 1D time vector, **x** and 2D data matrix, __y__.
+The function <code>plot_multichan()</code> draws multi-channel time series data, by taking 1D time vector, <code>x</code>, and 2D data matrix, <code>y</code>.
 
 ```python
 # Demo 2-2a. Function for multi-channel plotting
@@ -272,7 +275,8 @@ def plot_multichan( x, y, spacing = 3000, figsize = (10,10), ch_names = EEG.ch_n
 ```
   
   
-Using *plot_multichan()* function, example single-trial EEG trace can be visualized as follow.
+Using <code>plot_multichan()</code> function, example single-trial EEG trace can be visualized as follow.
+
 
 ```
 # Demo 2-2b. Visualization of raw EEG time trace
@@ -287,7 +291,7 @@ Note that channel 1 to 36 contain actual EEG data from 36-channel electrodes (fr
 
 ### 2-3. ERP in time domain
 
-Using *plot_multichan()*, ERP (Event-related potentials) trace can be drawn as follow.
+Using <code>plot_multichan()</code>, ERP (Event-related potentials) trace can be drawn as follow.
 
 ```python
 # Demo 2-3. Visualization of ERP time trace
@@ -302,45 +306,44 @@ plt.title('Condition #%d: %s (n=%d)'%(targetCondition, condNames[targetCondition
 
 ### 2-4. ERP in frequency domain
 
-To calculate the amplitude of 40-Hz auditory steady-state response, fast Fourier transform can be applied as follow. 
+As the sensory stimuli used in this experiment was 40 Hz click train, one can easily see the 40 Hz evoked (steady-state) response in the time-domain ERP. Through 
 
 ```python
-# Demo2-4.Time- and frequency-domain visualization of ERP
-fromnumpy.fftimportfft
-deffft_half(x,Fs=2000):returnfft(x)[:int(len(x)/2)]/len(x),np.linspace(0,Fs/2,len(x)/2)
+# Demo 2-4. Time- and frequency-domain visualization of grand-averaged ERP
+from numpy.fft import fft
+def fft_half(x, Fs=2000): return fft(x)[:int(len(x)/2)]/len(x), np.linspace(0,Fs/2,len(x)/2)
 
-# Get ERP in frequency domain
-erp=np.nanmean(EEG.data[:,:,np.where((EEG.events[:,2])==2)[0]],2)
-ac_erp=np.mean(erp[:36,],0)
-ac_erp_fft,freq=fft_half(ac_erp)
-
-# Plot time-domain ERP
 plt.figure(figsize=(14,8))
+
+erp = np.nanmean(EEG.data[:,:,np.where((EEG.events[:,2])==2)[0]],2)
+ac_erp = np.mean(erp[:36,],0)
+ac_erp_fft,freq = fft_half(ac_erp)
+
+# Plot time-domain signal
 plt.subplot(2,1,1)
-plt.plot(EEG.times,ac_erp,linewidth=.5)
+plt.plot( EEG.times, ac_erp, linewidth=.5 )
 plt.xlim((-.8,1.8))
-plt.xlabel(‘Time(sec)’)
-plt.ylabel(‘Amplitude(mV)’)
-plt.title(‘ERPsignalintimedomain’)
+plt.xlabel('Time (sec)')
+plt.ylabel('Amplitude (mV)')
+plt.title('ERP signal in time domain')
 plt.gca().set_facecolor((1,1,1))
 
-# Plot frequency-domain ERP
+# Plot frequency-domain signal
 plt.subplot(2,1,2)
-plt.plot(freq,np.abs(ac_erp_fft),linewidth=1)
+plt.plot( freq, np.abs(ac_erp_fft), linewidth=1 )
 plt.xlim((0,80))
-plt.xlabel(‘Freq(Hz)’)
-plt.ylabel(‘Amplitude(mV/Hz)’)
-plt.title(‘ERPsignalinfrequencydomain’)
+plt.xlabel('Freq (Hz)')
+plt.ylabel('Amplitude (mV/Hz)')
+plt.title('ERP signal in frequency domain')
 plt.gca().set_facecolor((1,1,1))
 
-plt.subplots_adjust(wspace=.3,hspace=.5)
-
+plt.subplots_adjust(wspace=.3, hspace=.5)
 ```
 
 
 ### 2-5. ERP in time-frequency domain
 
-Applying fast Fourier transform with moving temporal window, ERP signal can be drawn in time-frequency domain. To calculate spectrogram, a function *get_spectrogram()* is defined.
+Applying fast Fourier transform with moving temporal window, ERP signal can be drawn in time-frequency domain. To calculate spectrogram, a function <code>get_spectrogram()</code> is defined.
 
 ```python
 # Demo 2-5. Visualize frequency components in ERP
@@ -428,9 +431,9 @@ plt.subplots_adjust(wspace=.3, hspace=.5)
 <br>
 ### 3-1. Drawing 2D power topography
 
-Using this coordinate, spatial dynamics of ERP can be drawn in 2D plane (i.e., power topography). For visualization purpose, an example function, plot_topo2d( x, y ), is provided which takes x as 1D data matrix (i.e., 1 x 36 channels) each of which represents the channel power and takes y as 2D coordinate matrix. 
+Using this coordinate, spatial dynamics of ERP can be drawn in 2D plane (i.e., power topography). For visualization purpose, an example function, <code>plot_topo2d( data )</code>, is provided which takes <code>data</code> as 1D data matrix (i.e., 1 x 36 channels) each of which represents the channel power. 
 
-For 2D interpolation, additional *class* of *bi_interp2()* is defined. 
+For 2D interpolation, additional *class* of <code>bi_interp2</code> is defined. 
 
 ```python
 # Demo 3-1a. Preparation of 2D power topography
@@ -573,11 +576,12 @@ def plot_topo2d(data, clim=(-15,25), montage_file='%s%smontage.csv'%(dir_origin,
 ```
 
 **Optional) Bad channel identification**
- 
-In typical EEG recordings, large amplitude signals coming from bad channel usually be a problem. To prevent this, researchers usually eliminate bad channel data before drawing power topography. Here, simple bad channel identification method is implemented using channel correlation. Data from the bad channels identified here will be ignored and replaced by median value hereafter.
+
+In typical EEG recordings, large amplitude artifacts coming from few bad channel(s) usually be a problem. To prevent this, researchers have developed various methods of artifact-rejection and bad-channel selection. Here, simple bad channel identification method is implemented using channel correlation. Data from the bad channels identified here will be ignored in topography and replaced by median value hereafter.
 
 ```python
 # Demo 3-1b. Identification of bad-channel using correlation
+""" Bad-channel selection (correlation based) """
 from scipy.stats import ttest_1samp as ttest
 
 ga_erp = np.nanmean(EEG.data[:36,:,:],2)
@@ -606,7 +610,7 @@ print('\nLow-correlated (bad) channels: %s'%(bad_channels))
 
 ### 3-2. Time-course of raw voltage topography
 
-Input data of EEG topography can be defined by any mean; voltage, band-limited power, instantaneous angle, and so on. In this example, spatial distribution of raw voltage at specific time point is drawn. For better understanding of the data, ERP time traces at frontal and parietal area are also drawn.
+Input data of EEG topography can be defined by any mean; voltage, band-limited power, instantaneous phase, and so on. In this example, spatial distribution of raw voltage at specific time point is drawn. For better understanding of the data, ERP time traces at frontal and parietal area are also drawn.
 
 ```python
 # Demo 3-2. Raw voltage topography
@@ -624,29 +628,29 @@ for targetCondition in conditions:
   plt.subplot(len(conditions),1,np.where(np.array(conditions)==targetCondition)[0]+1)
   color_f = (.68,.210,.27) # Custom color value
   color_p = (.01,.457,.74)
-  plt.grid(‘off’)
-  plt.plot((0,0),(-30,45), ‘—‘, linewidth = 1, color = (.8,.8,.8))
-  plt.plot((-.2,1), (0,0), ‘—‘, linewidth = 1, color = (.8,.8,.8))
+  plt.grid('off')
+  plt.plot((0,0),(-30,45), '--', linewidth = 1, color = (.8,.8,.8))
+  plt.plot((-.2,1), (0,0), '--', linewidth = 1, color = (.8,.8,.8))
   plt.plot( EEG.times, frontal_erp,  color= color_f)
   plt.plot( EEG.times, parietal_erp, color= color_p)
-  plt.xlabel(‘Time (msec)’)
+  plt.xlabel('Time (msec)')
   plt.xlim((-.2,1))
   plt.ylim((-30,45))
-  plt.axis(‘off’)
+  plt.axis('off')
   plt.gca().set_facecolor((1,1,1))
-  plt.text( -.1, 27, ‘Frontal’, ha=‘center’, weight=‘bold’, fontsize=12, color=color_f )
-  plt.text( -.1, 20,’Parietal’, ha=‘center’, weight=‘bold’, fontsize=12, color=color_p )
-  plt.title(‘%s’%condNames[targetCondition-1])
+  plt.text( -.1, 27, 'Frontal', ha='center', weight='bold', fontsize=12, color=color_f )
+  plt.text( -.1, 20,'Parietal', ha='center', weight='bold', fontsize=12, color=color_p )
+  plt.title('%s'%condNames[targetCondition-1])
   
   # Calculate topography data 
   t_slice = [  (.005, .025), (.0250, .0550), (.200, .210),(.502, .512), (.925, .935)  ]  
   y_mark = [-20, 33, 10, 5, 10]
-  colors = [‘r’,’g’,’b’,’c’,’m’] # color marker
+  colors = ['r','g','b','c','m'] # color marker
   topos = []
   for tIdx in range(len(t_slice)):
     x_start, x_end, y_pos = t_slice[tIdx][0], t_slice[tIdx][1], y_mark[tIdx]
     idx_start, idx_end = np.where(EEG.times==x_start)[0][0], np.where(EEG.times==x_end)[0][0]
-    plt.plot( EEG.times[[idx_start,idx_end]], [y_pos, y_pos], colors[tIdx]+’|-‘)
+    plt.plot( EEG.times[[idx_start,idx_end]], [y_pos, y_pos], colors[tIdx]+'|-')
     topo_in = np.mean( erp[:36,idx_start:idx_end],1 )
     # bad-channel replacement
     topo_in[bad_channels] = np.median( topo_in.flatten() ) 
@@ -669,7 +673,7 @@ plt.close()
 
 ### 3-3. Band-limited power topography
 
-Other than raw voltage, topography of band-limited power at stimulation frequency (40 Hz) can be drawn as well. In this example, stimulus-evoked 40 Hz power were estimated using bandpower() function.  
+Other than raw voltage, topography of band-limited power at stimulation frequency (40 Hz) can be drawn as well. In this example, stimulus-evoked 40 Hz power were estimated using <code>bandpower()</code> function.  
 
 To demonstrate the effect of stimulation, stimulus-free periods (e.g., pre- and post-stimulus period) data are also obtained.
 
@@ -692,20 +696,20 @@ targetCondition = 6 # = Auditory sound only
 trialIdx = np.where((EEG.events[:,2])==targetCondition)[0]
 erp = np.mean(EEG.data[:,:,trialIdx],2)
 period = [ (-.5,0.), (0.,.5), (.5, 1.), (1.,1.5) ] # time in second
-periodName = [‘Pre-stim’, ‘Stim (Early)’, ‘Stim (Late)’, ‘Post-stim’];
+periodName = ['Pre-stim', 'Stim (Early)', 'Stim (Late)', 'Post-stim'];
 freq = 40 # Hz
 plt.figure(figsize=(15,3))
 for periodIdx in range(len(period)):
   tIdx = (EEG.times>period[periodIdx][0]) & (EEG.times<=period[periodIdx][1])
   
   # Calculate power & Substitute bad-channel value
-  power = band_power(erp[:36,tIdx], np.array([-2,2])+freq, EEG.info[‘sfreq’])
+  power = band_power(erp[:36,tIdx], np.array([-2,2])+freq, EEG.info['sfreq'])
   power[bad_channels]= np.median(power.flatten())
   
   # Draw
   plt.subplot(1,len(period),periodIdx+1)
   plot_topo2d(power, clim=(0,3) )
-  plt.title(‘%s, t = [%.1f, %.1f]s’%(periodName[periodIdx],period[periodIdx][0],period[periodIdx][1]))
+  plt.title('%s, t = [%.1f, %.1f]s'%(periodName[periodIdx],period[periodIdx][0],period[periodIdx][1]))
 
 ```
 
@@ -725,15 +729,13 @@ for targetCondition in conditions:
   erp = np.mean(EEG.data[:,:,trialIdx],2)
   
   # Calculate power & Substitute bad-channel value
-  power = band_power(erp[:36,tIdx], np.array([-2,2])+freq, EEG.info[‘sfreq’])
+  power = band_power(erp[:36,tIdx], np.array([-2,2])+freq, EEG.info['sfreq'])
   power[bad_channels]= np.median(power.flatten())
   
   # Draw
   plt.subplot(1,len(conditions),np.where(np.array(conditions)==targetCondition)[0]+1)
   plot_topo2d(power, clim=(0,7) )
-  plt.title(‘%s’%condNames[targetCondition-1])
-  if targetCondition is not conditions[0]: plt.ylabel(‘’)
+  plt.title('%s'%condNames[targetCondition-1])
+  if targetCondition is not conditions[0]: plt.ylabel('')
 ```
-<br>
 
-Enjoy!