- "subs = [101, 102, 103, 104, 105, 106, 108, 109, 110, 111, 112,\n 202, 203, 204, 205, 207, 208, 209, 210, 211, \n 301, 302, 303, 304, 305, 306, 307, 308, 309]\n\ndiff_out = []\nIpsi_out = []\nContra_out = []\nIpsi_spectra_out = []\nContra_spectra_out = []\ndiff_spectra_out = []\nERSP_diff_out = []\nERSP_Ipsi_out = []\nERSP_Contra_out = []\n\nfrequencies = np.linspace(6, 30, 100, endpoint=True)\nwave_cycles = 6\n\n# time frequency window for analysis\nf_low = 7 # Hz\nf_high = 10\nf_diff = f_high-f_low\n \nt_low = 0 # s\nt_high = 1\nt_diff = t_high-t_low\n\nbad_subs= [6, 7, 13, 26]\nreally_bad_subs = [11, 12, 19]\nsub_count = 0 \n \n \n \nfor sub in subs:\n print(sub)\n \n sub_count += 1\n\n \n if (sub_count in really_bad_subs):\n rej_thresh_uV = 90\n elif (sub_count in bad_subs):\n rej_thresh_uV = 90\n else:\n rej_thresh_uV = 90\n\n rej_thresh = rej_thresh_uV*1e-6\n \n \n # Load both sessions\n raw = load_data(sub,1, # subject, session\n experiment='visual-cueing',site='kylemathlab_dev',device_name='muse2016',\n data_dir = eegnb_data_path)\n \n raw.append(\n load_data(sub,2, # subject, session\n experiment='visual-cueing', site='kylemathlab_dev', device_name='muse2016',\n data_dir = eegnb_data_path))\n \n\n # Filter Raw Data\n raw.filter(1,30, method='iir')\n\n #Select Events\n events = find_events(raw)\n event_id = {'LeftCue': 1, 'RightCue': 2}\n epochs = Epochs(raw, events=events, event_id=event_id, \n tmin=-1, tmax=2, baseline=(-1, 0), \n reject={'eeg':rej_thresh}, preload=True,\n verbose=False, picks=[0, 3])\n print('Trials Remaining: ' + str(len(epochs.events)) + '.')\n\n # Compute morlet wavelet\n\n # Left Cue\n tfr, itc = tfr_morlet(epochs['LeftCue'], freqs=frequencies, \n n_cycles=wave_cycles, return_itc=True)\n tfr = tfr.apply_baseline([-1,-.5],mode='mean')\n #tfr.plot(picks=[0], mode='logratio', \n # title='TP9 - Ipsi');\n #tfr.plot(picks=[3], mode='logratio', \n # title='TP10 - Contra');\n power_Ipsi_TP9 = tfr.data[0,:,:]\n power_Contra_TP10 = tfr.data[1,:,:]\n\n # Right Cue\n tfr, itc = tfr_morlet(epochs['RightCue'], freqs=frequencies, \n n_cycles=wave_cycles, return_itc=True)\n tfr = tfr.apply_baseline([-1,-.5],mode='mean')\n #tfr.plot(picks=[0], mode='logratio', \n # title='TP9 - Contra');\n #tfr.plot(picks=[3], mode='logratio', \n # title='TP10 - Ipsi');\n power_Contra_TP9 = tfr.data[0,:,:]\n power_Ipsi_TP10 = tfr.data[1,:,:]\n\n # Plot Differences\n #%matplotlib inline\n times = epochs.times\n power_Avg_Ipsi = (power_Ipsi_TP9+power_Ipsi_TP10)/2;\n power_Avg_Contra = (power_Contra_TP9+power_Contra_TP10)/2;\n power_Avg_Diff = power_Avg_Ipsi-power_Avg_Contra;\n\n\n #find max to make color range\n plot_max = np.max([np.max(np.abs(power_Avg_Ipsi)), np.max(np.abs(power_Avg_Contra))])\n plot_diff_max = np.max(np.abs(power_Avg_Diff))\n\n \n \n #Ipsi\n fig, ax = plt.subplots(1)\n im = plt.imshow(power_Avg_Ipsi,\n extent=[times[0], times[-1], frequencies[0], frequencies[-1]],\n aspect='auto', origin='lower', cmap='coolwarm', vmin=-plot_max, vmax=plot_max)\n plt.xlabel('Time (sec)')\n plt.ylabel('Frequency (Hz)')\n plt.title('Power Average Ipsilateral to Cue')\n cb = fig.colorbar(im)\n cb.set_label('Power')\n # Create a Rectangle patch\n rect = patches.Rectangle((t_low,f_low),t_diff,f_diff,linewidth=1,edgecolor='k',facecolor='none')\n # Add the patch to the Axes\n ax.add_patch(rect)\n\n #TP10\n fig, ax = plt.subplots(1)\n im = plt.imshow(power_Avg_Contra,\n extent=[times[0], times[-1], frequencies[0], frequencies[-1]],\n aspect='auto', origin='lower', cmap='coolwarm', vmin=-plot_max, vmax=plot_max)\n plt.xlabel('Time (sec)')\n plt.ylabel('Frequency (Hz)')\n plt.title(str(sub) + ' - Power Average Contra to Cue')\n cb = fig.colorbar(im)\n cb.set_label('Power')\n # Create a Rectangle patch\n rect = patches.Rectangle((t_low,f_low),t_diff,f_diff,linewidth=1,edgecolor='k',facecolor='none')\n # Add the patch to the Axes\n ax.add_patch(rect)\n\n #difference between conditions\n fig, ax = plt.subplots(1)\n im = plt.imshow(power_Avg_Diff,\n extent=[times[0], times[-1], frequencies[0], frequencies[-1]],\n aspect='auto', origin='lower', cmap='coolwarm', vmin=-plot_diff_max, vmax=plot_diff_max)\n plt.xlabel('Time (sec)')\n plt.ylabel('Frequency (Hz)')\n plt.title('Power Difference Ipsi-Contra')\n cb = fig.colorbar(im)\n cb.set_label('Ipsi-Contra Power')\n # Create a Rectangle patch\n rect = patches.Rectangle((t_low,f_low),t_diff,f_diff,linewidth=1,edgecolor='k',facecolor='none')\n # Add the patch to the Axes\n ax.add_patch(rect)\n \n \n \n \n #output data into array\n Ipsi_out.append(np.mean(power_Avg_Ipsi[np.argmax(frequencies>f_low):\n np.argmax(frequencies>f_high)-1,\n np.argmax(times>t_low):np.argmax(times>t_high)-1 ]\n )\n ) \n Ipsi_spectra_out.append(np.mean(power_Avg_Ipsi[:,np.argmax(times>t_low):\n np.argmax(times>t_high)-1 ],1\n )\n )\n \n Contra_out.append(np.mean(power_Avg_Contra[np.argmax(frequencies>f_low):\n np.argmax(frequencies>f_high)-1,\n np.argmax(times>t_low):np.argmax(times>t_high)-1 ]\n )\n )\n \n Contra_spectra_out.append(np.mean(power_Avg_Contra[:,np.argmax(times>t_low):\n np.argmax(times>t_high)-1 ],1))\n \n \n diff_out.append(np.mean(power_Avg_Diff[np.argmax(frequencies>f_low):\n np.argmax(frequencies>f_high)-1,\n np.argmax(times>t_low):np.argmax(times>t_high)-1 ]\n )\n )\n diff_spectra_out.append(np.mean(power_Avg_Diff[:,np.argmax(times>t_low):\n np.argmax(times>t_high)-1 ],1\n )\n )\n \n \n ERSP_diff_out.append(power_Avg_Diff)\n ERSP_Ipsi_out.append(power_Avg_Ipsi)\n ERSP_Contra_out.append(power_Avg_Contra)\n\n\n \nprint(np.shape(ERSP_diff_out))\nprint(np.shape(Contra_spectra_out))\n\nprint(diff_out)"
0 commit comments