OSCR

Fixation-related potentials reveal that confusing program code elicits a late frontal positivity.

Code ↔ Paper

5 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 5 matches · 1 of them tie a paragraph to a whole file, not to given lines: a weak match, whose lines are not tinted
  1. [1] § Methods › Recording and preprocessing › EEG preprocessing ↔ 09-Task-Evaluate_Data/utils/eeg_helpers.py, lines 572–686 · score 0.72 · Brain Vision, EEG channels, weight, segments, algorithm
  2. [2] § Methods › Recording and preprocessing › Eye-movement recording ↔ 06-Task-Study_Presentation_Software/Experiment_lastrun.py, lines 81–97 · score 0.71 · Tobii Pro Spectrum, Tobii Eye Tracker
  3. [3] § Methods › Recording and preprocessing › Eye-movement preprocessing ↔ 09-Task-Evaluate_Data/utils/I2MC_settings.py, the whole file · a weak match · score 0.70 · I2MC, fixation duration, eye tracking, noise, saccades, algorithm
  4. [4] § Methods › Data analysis › FRP analysis ↔ 09-Task-Evaluate_Data/utils/eeg_helpers.py, lines 755–829 · score 0.67 · baseline corrected, stimulus onset, voltage, segment, amplitude, absolute
  5. [5] § Methods › Recording and preprocessing › EEG preprocessing ↔ 09-Task-Evaluate_Data/utils/path_helpers.py, lines 309–317 · score 0.64 · Brain Vision Analyzer, EEG

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

Python · 1,133 lines · 65 KB · CC-BY-4.0 · 2 matches

  1. import gc
  2. import json
  3. import re
  4. from math import log10
  5. from pathlib import Path
  6. from typing import Union
  8. import mne
  9. import numpy as np
  10. import pandas as pd
  11. import seaborn as sns
  12. from matplotlib import pyplot as plt
  13. import matplotlib as mpl
  14. from mne.io import Raw
  15. from tqdm.notebook import tqdm
  16. from utils.eeg_settings import (ACCEPTED_SYNCHRONIZATION_OFFSET, EEG_CHANNELS,
  17. EEG_FREQUENCY, EEG_LONG_BUFFER,
  18. EEG_MEAN_BUFFER, EEG_SHORT_BUFFER,
  19. EEG_STIMULUS, EEG_STIMULUS_FIXATION_CROSS,
  20. EEG_STIMULUS_SNIPPET_END,
  21. EEG_STIMULUS_SNIPPET_START,
  22. EEG_VOLTAGE_OVERALL, EEG_VOLTAGE_STEP,
  23. EEG_VOLTAGE_WINDOW, EOG_CHANNELS, ERP_PARAMETER_CORRECT_TRIALS_ONLY, ERP_PARAMETER_EPOCH_INTERVAL,
  24. FRP_EEG_STIMULUS_SNIPPET_START,
  25. IMPEDANCE_UPPER_BOUND, IMPEDANCE_VALUE, MNE_KEY_FREQUENCY,
  26. STIMULUS_EVENT_NAMES)
  27. from utils.file_helpers import (get_exclusions,
  28. get_participant_folder_per_participant)
  29. from utils.file_settings import (ANNOTATION_COLUMN_DESCRIPTION,
  30. ANNOTATION_COLUMN_ONSET,
  31. ANNOTATION_COLUMN_ONSET_FLOAT,
  32. BEHAVIORAL_COLUMN_CORRECTNESS,
  33. BEHAVIORAL_COLUMN_END,
  34. BEHAVIORAL_COLUMN_FIXATION_START,
  35. BEHAVIORAL_COLUMN_START, COLUMN_TIME,
  36. EEG_COLUMN_STIMULUS, FIXATION_COLUMN_START,
  37. HDF_INDEX, SEPARATOR)
  38. from utils.path_helpers import (get_all_erp_epoch_paths, get_behavioral_data_path, get_erp_average_path,
  39. get_eeg_trial_path, get_erp_epoch_path,
  40. get_erp_fixation_analysis_path,
  41. get_erp_nave_path, get_erp_status_path)
  42. from utils.path_settings import (EEG_FILE_DATA_ENDING, EEG_FILE_HEADER_ENDING,
  43. EEG_FILE_MARKER_ENDING, PROCESSED_PATH)
  44. from utils.snippet_helpers import get_snippet_number, get_snippet_variant
  45. from utils.snippet_settings import (CONDITION, CONDITION_CLEAN,
  46. CONDITION_COLORS, CONDITION_CONFUSING,
  47. CONDITION_DIFF, CONDITION_VARIANT_MATCH,
  48. PANDAS_DESCRIPTION_AGG_FUNCTIONS,
  49. PANDAS_DESCRIPTION_AGG_NAMES, SNIPPET_GROUP_ALL, SNIPPET_NUMBERS)
  50. from utils.textconstants import (BEHAVIORAL, EEG, EEG_ERP, FIXATIONS,
  51. PARTICIPANT, SNIPPET, TIME, TOTAL, VISUAL)
  52. from utils.visual_settings import (FIXATION_SELECTION_ALGORITHM,
  53. FIXATION_SELECTION_ALGORITHMS,
  54. FIXATION_SELECTION_SHORT_VERSION)
  57. def check_file_existence(files: dict[str, Path], file: Path, file_ending: str, participant: str):
  58. '''check whether file has the given extension and there already exists one.
  60. Arguments:
  61. * files: where to add suitable files per ending
  62. * file: the path to check
  63. * file_endings: the file ending to check for
  64. * participant: the participant to name when problems arise
  66. raises: Exception if file with given ending has already been identified
  67. '''
  68. if file.suffix == file_ending:
  69. if file_ending in files:
  70. print(
  71. f'Multiple \'{file_ending}\' files found for participant {participant}: {files[file_ending].name}, new: {file.name}')
  72. raise Exception()
  73. files[file_ending] = file
  76. def get_eeg_files_per_participant() -> dict[str, dict[str, Path]]:
  77. f'''get eeg files (3 files as dictionary) per participant folders identified in the base path
  78. requirement: exactly one eeg file of each file ending per participant (or none at all, then participant is ignored)
  80. returns: the participant numbers, and per each the three eeg files
  81. {{'{EEG_FILE_HEADER_ENDING}': path to the eeg header file,
  82. '{EEG_FILE_MARKER_ENDING}': path to the eeg marker file,
  83. '{EEG_FILE_DATA_ENDING}': path to the eeg data file}}
  84. '''
  85. eeg_files = {}
  87. for participant, participant_folder in get_participant_folder_per_participant().items():
  88. files = {}
  89. for file in participant_folder.iterdir():
  90. check_file_existence(
  91. files, file, EEG_FILE_DATA_ENDING, participant)
  92. check_file_existence(
  93. files, file, EEG_FILE_HEADER_ENDING, participant)
  94. check_file_existence(
  95. files, file, EEG_FILE_MARKER_ENDING, participant)
  96. # check that either all three files exits, or none at all
  97. if ((EEG_FILE_DATA_ENDING in files) ^ (EEG_FILE_HEADER_ENDING in files)) or ((EEG_FILE_DATA_ENDING in files) ^ (EEG_FILE_MARKER_ENDING in files)):
  98. print(f'Not all files found for {participant}: {files}')
  99. return None
  100. # only add if all files exist
  101. if files:
  102. eeg_files[participant] = files
  103. return eeg_files
  106. EEG_HEADER_DATE = re.compile('Impedance \[kOhm\] at (\d\d:\d\d:\d\d) :')
  107. EEG_MARKER_DATE = re.compile(r'New Segment,,(\d+),1,0,(\d+)')
  110. def anonymize_eeg_data(eeg_file: dict[str, Path]):
  111. '''anonymizes eeg files in-place by removing all traces of timestamps in the marker and header files
  113. Argument: eeg_file: a dictionary mapping the eeg file keys to the respective paths of the files
  114. '''
  115. with open(eeg_file[EEG_FILE_HEADER_ENDING], 'r+') as f:
  116. eeg_header_content = f.read()
  117. f.seek(0)
  118. eeg_header_content = EEG_HEADER_DATE.sub(
  119. 'Impedance [kOhm] at the beginning of the experiment :', eeg_header_content)
  120. f.write(eeg_header_content)
  121. f.truncate()
  122. with open(eeg_file[EEG_FILE_MARKER_ENDING], 'r+') as f:
  123. eeg_marker_content = f.read()
  124. f.seek(0)
  125. eeg_marker_content = EEG_MARKER_DATE.sub(
  126. r'New Segment,,\1,1,0,0', eeg_marker_content, count=0)
  127. f.write(eeg_marker_content)
  128. f.truncate()
  131. def load_eeg_data(filepath: Path, preload: bool = True) -> tuple[Raw, float]:
  132. raw_eeg = mne.io.read_raw_brainvision(filepath, eog=tuple(EOG_CHANNELS),
  133. preload=preload)
  134. frequency = raw_eeg.info[MNE_KEY_FREQUENCY]
  135. return raw_eeg, frequency
  138. def check_impedances(impedance_data: pd.DataFrame, log: bool = True) -> bool:
  139. impedance_data[f'{IMPEDANCE_UPPER_BOUND}_check'] = impedance_data.apply(lambda row: max(
  140. 0, row[IMPEDANCE_VALUE]-row[IMPEDANCE_UPPER_BOUND]), 1)
  141. # display(impedance_data)
  142. impedance_errors = ((
  143. impedance_data[f'{IMPEDANCE_UPPER_BOUND}_check'] > 0)*1).sum(), \
  144. impedance_data[f'{IMPEDANCE_UPPER_BOUND}_check'].max(), \
  145. impedance_data[f'{IMPEDANCE_UPPER_BOUND}_check'].sum()
  147. impedance_okay = not (impedance_errors[0] > 2 or
  148. impedance_errors[1] > 5 or
  149. impedance_errors[2] >= 7)
  150. if log:
  151. print(
  152. f'{impedance_errors[0]} errors with impedances, maximum breach of {impedance_errors[1]}, sum of breaches in total {impedance_errors[2]}.\n\tImpedances accepted: {impedance_okay}')
  153. return impedance_okay
  156. def prepare_annotation_information(eeg_data: Raw) -> tuple[pd.DataFrame, float]:
  157. '''extracts and prepares annotation information from the given EEG data
  159. Arguments: eeg_data: the eeg data
  161. returns: a tuple of
  162. * the annotations as DataFrame with columns onset, duration and description extracted from the eeg data, as well as an additional column describing the onset as a float value to be used for cropping, reinsert, ...
  163. * the duration of the eeg data in seconds
  164. '''
  165. # get recording information required+
  166. # offset in seconds between start of the file counter and start of samples
  167. recording_offset = eeg_data.first_time
  168. recording_duration = (eeg_data.n_times-1) / \
  169. eeg_data.info[MNE_KEY_FREQUENCY] # duration of recording
  171. annotation_data = eeg_data.annotations.to_data_frame()
  172. # transform onset of an annotation into a float based on recording start, erasing offset, to use for cropping
  173. annotation_data[ANNOTATION_COLUMN_ONSET_FLOAT] = annotation_data[ANNOTATION_COLUMN_ONSET].apply(
  174. lambda onset: (onset-pd.Timestamp(year=1970, month=1, day=1)).total_seconds()-recording_offset)
  176. return annotation_data, recording_duration
  179. def crop_to_complete_annotation_range(eeg_data: Raw) -> None:
  180. '''crop given Raw object by the annotations present in the object.
  182. The new sequence starts shortly before first fixation cross,
  183. and end shortly after last snippet end, or longer after the last fixation cross
  185. Arguments:
  186. * eeg_data: the eeg data
  187. * time_before_first_snippet: the time buffer to add before the first fixation cross
  188. * time_after_last_snippet: the time buffer to add after the last snippet end if found
  189. * time_constant_without_ending: the time buffer to add after the last fixation cross if no end found
  190. '''
  191. # get required information from annotations
  192. annotation_data, recording_duration = prepare_annotation_information(
  193. eeg_data)
  195. # get first fixation cross or start
  196. start_buffer = {EEG_STIMULUS_FIXATION_CROSS: EEG_SHORT_BUFFER} | {
  197. stimuli: EEG_MEAN_BUFFER for stimuli in EEG_STIMULUS_SNIPPET_START.values()}
  198. snippets_start_row = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin(
  199. start_buffer.keys())].iloc[0]
  201. snippets_start_time = snippets_start_row[ANNOTATION_COLUMN_ONSET_FLOAT] - \
  202. start_buffer[snippets_start_row[ANNOTATION_COLUMN_DESCRIPTION]]
  203. snippets_start_time = max(0.0, snippets_start_time)
  205. # get snippet ends and use them to calculate the end of the cropped recording
  206. end_buffer = {EEG_STIMULUS_SNIPPET_END: EEG_SHORT_BUFFER} | {
  207. stimuli: EEG_LONG_BUFFER for stimuli in EEG_STIMULUS_SNIPPET_START.values()}
  208. snippets_end_row = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin([
  209. EEG_STIMULUS_FIXATION_CROSS, *EEG_STIMULUS_SNIPPET_START.values(), EEG_STIMULUS_SNIPPET_END])]
  211. if not snippets_end_row.empty:
  212. snippets_end_time = snippets_end_row.iloc[-1][ANNOTATION_COLUMN_ONSET_FLOAT] + \
  213. EEG_MEAN_BUFFER
  214. # otherwise, add a longer buffer to the last fixation
  215. else:
  216. snippets_end_time = min(
  217. recording_duration, snippets_end_row[ANNOTATION_COLUMN_ONSET_FLOAT] +
  218. end_buffer[snippets_end_row[ANNOTATION_COLUMN_DESCRIPTION]])
  219. # crop
  220. assert (recording_duration >= snippets_end_time-snippets_start_time)
  221. eeg_data.crop(tmin=snippets_start_time,
  222. tmax=snippets_end_time, include_tmax=True)
  225. def export_eeg_brainvision(eeg_data: Raw, eeg_path: Path):
  226. assert (eeg_path.suffix == EEG_FILE_HEADER_ENDING)
  227. eeg_data.export(eeg_path, overwrite=True, verbose=False)
  228. eeg_marker_path = eeg_path.with_suffix(EEG_FILE_MARKER_ENDING)
  229. with open(eeg_marker_path, 'r+') as f:
  230. eeg_marker_content = f.read()
  231. f.seek(0)
  232. eeg_marker_content = eeg_marker_content.replace(
  233. 'Comment,Bad Interval/', 'Bad Interval,')
  234. eeg_marker_content = eeg_marker_content.replace(
  235. 'Comment,Bad', 'Bad Interval,Bad')
  236. f.write(eeg_marker_content)
  237. f.truncate()
  240. def check_manual_ICA_reasoning(artifact_reasoning: pd.DataFrame):
  241. # print(artifact_reasoning)
  242. # check that all components are still there
  243. expected_components = set(f'F{str(i).zfill(2)}' for i in range(32))
  244. given_components = set(artifact_reasoning['Component'].values)
  245. assert given_components == expected_components, f'''The components are not correct, \n\twanted {
  246. expected_components},\n\tgiven {given_components}'''
  247. # check that description & topology are filled
  248. assert (artifact_reasoning['Description'].apply(lambda v: not pd.isna(
  249. v)).all()), 'Description must be filled for each component.'
  250. assert (artifact_reasoning['Topology'].apply(lambda v: not pd.isna(
  251. v)).all()), 'Topology must be filled for each component.'
  252. # check that reason is given if artifact... is not false
  253. assert artifact_reasoning.apply(lambda row: not pd.isna(row['Reason']) if row['Artifact or Channel related'] != False else True, axis=1).all(
  254. ), 'Each component identified as possibly being artifact or channel related must have a reason towards the choice of inclusion or not.'
  255. # check that included is not false if artifact... is false
  256. assert artifact_reasoning.apply(lambda row: row['Included'] != False if row['Artifact or Channel related'] == False else True, axis=1).all(
  257. ), 'Each component not identified as possibly being artifact or channel related must be included.'
  260. def assign_trials_to_annotations(eeg_data: Raw, behavioral_events: pd.DataFrame) -> tuple[bool, pd.DataFrame]:
  261. '''Split given Raw object to receive eeq splits per snippet.
  263. Each split starts shortly before first fixation cross,
  264. and end shortly after next snippet end, or until the next fixation cross.
  265. If there is no other way of determination, a long buffer is added to the current fixation crops of the split.
  267. Arguments:
  268. * eeg_data: the eeg data
  269. * sequence_order: the snippet sequence order to assign eeg splits to each snippet
  270. * time_before_first_snippet: the time buffer to add before the first fixation cross
  271. * time_after_last_snippet: the time buffer to add after the last snippet end if found
  272. * time_constant_without_ending: the time buffer to add after the last fixation cross if no end found
  274. returns: eeg split per snippet
  275. '''
  276. annotation_data, _ = prepare_annotation_information(
  277. eeg_data)
  278. annotation_data = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin([EEG_STIMULUS_FIXATION_CROSS,
  279. EEG_STIMULUS_SNIPPET_END,
  280. *EEG_STIMULUS_SNIPPET_START.values()])]
  282. def annotation_synchronization_check(data: pd.DataFrame, time_check: bool = False) -> bool:
  283. data['stimuli_check'] = data[ANNOTATION_COLUMN_DESCRIPTION] == data[EEG_COLUMN_STIMULUS]
  284. if time_check:
  285. data['time_check'] = data.apply(lambda row: (row['onset_e'] < 0) or
  286. ((row['onset_a'] - row['onset_e']) < ACCEPTED_SYNCHRONIZATION_OFFSET), axis=1)
  287. if not data['stimuli_check'].all():
  288. return False
  289. if time_check and not data['time_check'].all():
  290. return False
  291. data = data.drop(columns=[c for c in data if (c == 'onset_a') or (
  292. (c != SNIPPET) and (not c in annotation_data.columns))])
  293. return True
  295. print('annotations:', annotation_data.shape[0],
  296. 'behavioral_events:', behavioral_events.shape[0])
  298. annotation_data.reset_index(drop=True, inplace=True)
  299. behavioral_events.reset_index(drop=True, inplace=True)
  300. # possibility 1: match via index if both are the same length
  301. if annotation_data.shape[0] == behavioral_events.shape[0]:
  302. annotation_data['onset_a'] = 0
  303. behavioral_events['onset_e'] = 0
  304. for hdf_index in behavioral_events[HDF_INDEX].unique():
  305. hdf_data = behavioral_events[behavioral_events[HDF_INDEX] == hdf_index]
  306. anno_data = annotation_data[behavioral_events[HDF_INDEX] == hdf_index]
  307. behavioral_events.loc[hdf_data.index,
  308. 'onset_e'] = hdf_data['Time'] - hdf_data['Time'].iloc[0]
  309. annotation_data.loc[anno_data.index, 'onset_a'] = anno_data[ANNOTATION_COLUMN_ONSET_FLOAT] - \
  310. anno_data[ANNOTATION_COLUMN_ONSET_FLOAT].iloc[0]
  312. com_data = pd.concat([behavioral_events, annotation_data], axis=1)
  313. if annotation_synchronization_check(com_data):
  314. return True, com_data
  315. test_data = pd.concat([behavioral_events, annotation_data], axis=1)
  316. # possibility 2: match via time based on first event / annotation
  317. # com_data = pd.DataFrame([], columns = [c for c in test_data.columns])
  318. # annotations_index, behavioral_index=0, 0
  319. # annotation_data['onset_a'] = 0
  320. # behavioral_events['onset_e'] = 0
  321. # for hdf_index in behavioral_events[HDF_INDEX].unique():
  322. # hdf_data = behavioral_events[behavioral_events[HDF_INDEX]==hdf_index]
  323. # anno_data = annotation_data[behavioral_events[HDF_INDEX]==hdf_index]
  324. # behavioral_events.loc[hdf_data.index, 'onset_e'] = hdf_data['Time'] - hdf_data['Time'].iloc[0]
  325. # annotation_data.loc[anno_data.index, 'onset_a'] = anno_data[ANNOTATION_COLUMN_ONSET_FLOAT] - anno_data[ANNOTATION_COLUMN_ONSET_FLOAT].iloc[0]
  326. # while (annotations_index<annotation_data.shape[0] and behavioral_index<behavioral_events.shape[0]):
  327. # pass
  328. # if (not com_data.empty) and annotation_synchronization_check(com_data):
  329. # return True, com_data
  330. # possibility 3: match via time based on last event / annotation
  331. # com_data = pd.DataFrame([], columns = [c for c in test_data.columns])
  332. # annotations_index, behavioral_index=annotation_data.shape[0]-1, behavioral_events.shape[0]-1
  333. # while (annotations_index>=0 and behavioral_index<behavioral_events>=0):
  334. # if (not com_data.empty) and annotation_synchronization_check(com_data):
  335. # return True, com_data
  336. # possibility 4: match via time diff based on events / annotation
  337. # com_data = pd.DataFrame([], columns = [c for c in test_data.columns])
  338. # if (not com_data.empty) and annotation_synchronization_check(com_data):
  339. # return True, com_data
  340. return False, test_data
  343. def check_trial_annotations(trial_annotation_data: pd.DataFrame, eeg_data: Raw):
  344. # print(trial_annotation_data)
  345. annotation_data, _ = prepare_annotation_information(eeg_data)
  346. annotation_data = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin([EEG_STIMULUS_FIXATION_CROSS,
  347. EEG_STIMULUS_SNIPPET_END,
  348. *EEG_STIMULUS_SNIPPET_START.values()])]
  349. data1 = (trial_annotation_data[[ANNOTATION_COLUMN_DESCRIPTION, ANNOTATION_COLUMN_ONSET_FLOAT]].dropna(axis='index')
  350. .sort_values([ANNOTATION_COLUMN_ONSET_FLOAT]).reset_index(drop=True))
  351. data2 = (annotation_data[[ANNOTATION_COLUMN_DESCRIPTION, ANNOTATION_COLUMN_ONSET_FLOAT]]
  352. .sort_values([ANNOTATION_COLUMN_ONSET_FLOAT]).reset_index(drop=True))
  353. # check that all annotations covered by a line
  354. assert data1[ANNOTATION_COLUMN_DESCRIPTION].equals(
  355. data2[ANNOTATION_COLUMN_DESCRIPTION]), 'All annotations must be present.'
  356. data1[ANNOTATION_COLUMN_ONSET_FLOAT] = data1[ANNOTATION_COLUMN_ONSET_FLOAT].round(
  357. 3)
  358. data2[ANNOTATION_COLUMN_ONSET_FLOAT] = data2[ANNOTATION_COLUMN_ONSET_FLOAT].round(
  359. 3)
  360. time_delta = (data2[ANNOTATION_COLUMN_ONSET_FLOAT] -
  361. data1[ANNOTATION_COLUMN_ONSET_FLOAT]).abs().ge(1.5/EEG_FREQUENCY)
  362. assert not time_delta.any(), f'All annotations must be present with their given frame'
  363. # check that all annotations, that are not fixation crosses, have an assigned snippet
  364. assert (trial_annotation_data.apply(lambda row: (pd.isna(row[ANNOTATION_COLUMN_DESCRIPTION]) or row[ANNOTATION_COLUMN_DESCRIPTION] == EEG_STIMULUS_FIXATION_CROSS) or pd.notna(
  365. row[SNIPPET]), axis=1).all()), "All annotations that are not fixation crosses require an assigned event."
  366. # check that stimuli identical
  367. assert (trial_annotation_data.dropna(axis='index').apply(lambda row: row[EEG_COLUMN_STIMULUS] == row[ANNOTATION_COLUMN_DESCRIPTION], axis=1).all(
  368. )), 'The stimuli of event and annotation must be identical'
  369. # check time synchronization between events and annotations work within a hdf file index
  370. for hdf_index in trial_annotation_data[HDF_INDEX].unique():
  371. if pd.isna(hdf_index):
  372. continue
  373. hdf_trials = trial_annotation_data[trial_annotation_data[HDF_INDEX] == hdf_index].dropna(
  374. axis='index')
  375. if hdf_trials.empty:
  376. continue
  377. hdf_trials['Time delta beh eeg'] = hdf_trials[ANNOTATION_COLUMN_ONSET_FLOAT] - hdf_trials[TIME]
  378. assert (hdf_trials['Time delta beh eeg'].max() - hdf_trials['Time delta beh eeg'].min() < 0.3), \
  379. f'''The difference between eeg annotation timestamp and behavioral event timestamp should remain within a second of time.\n {
  380. hdf_trials["Time delta beh eeg"]}'''
  383. def get_synchronized_annotations(trial_annotation_data: pd.DataFrame, behavioral_data: pd.DataFrame) -> pd.DataFrame:
  384. annotations_to_delete = []
  385. # check isna and print na lines to delete
  386. na_rows = trial_annotation_data[trial_annotation_data.isna().any(axis=1)]
  387. if not na_rows.empty:
  388. print(f'\tThese rows will be deleted (ignored) due to nas in the rows')
  389. print(na_rows)
  390. relevant_trial_annotation_data = trial_annotation_data.dropna(axis='index')
  391. for i, row in relevant_trial_annotation_data.iterrows():
  392. behavioral_row = behavioral_data[behavioral_data[SNIPPET] == row[SNIPPET]].squeeze(
  393. )
  394. if behavioral_row.empty:
  395. annotations_to_delete.append(i)
  396. continue
  397. if row[ANNOTATION_COLUMN_DESCRIPTION] == EEG_STIMULUS_FIXATION_CROSS:
  398. behavioral_time = behavioral_row[BEHAVIORAL_COLUMN_FIXATION_START]
  399. elif row[ANNOTATION_COLUMN_DESCRIPTION] in EEG_STIMULUS_SNIPPET_START.values():
  400. behavioral_time = behavioral_row[BEHAVIORAL_COLUMN_START]
  401. elif row[ANNOTATION_COLUMN_DESCRIPTION] == EEG_STIMULUS_SNIPPET_END:
  402. behavioral_time = behavioral_row[BEHAVIORAL_COLUMN_END]
  403. if abs(row[TIME] - behavioral_time) > 0.0001:
  404. annotations_to_delete.append(i)
  405. print(
  406. f'\tAnnotation {row} ignored even though assigned behavioral event, as the behavioral time is expected to be {behavioral_time} from {behavioral_row[[SNIPPET, PARTICIPANT, BEHAVIORAL_COLUMN_FIXATION_START, BEHAVIORAL_COLUMN_START, BEHAVIORAL_COLUMN_END]]}.')
  407. relevant_trial_annotation_data = relevant_trial_annotation_data[~relevant_trial_annotation_data.index.isin(
  408. annotations_to_delete)]
  409. return relevant_trial_annotation_data
  412. def transform_synchronized_annotations(trial_annotation_data: pd.DataFrame) -> pd.DataFrame:
  413. eeg_snippet_data = trial_annotation_data.pivot(
  414. columns=EEG_COLUMN_STIMULUS, index=SNIPPET, values=ANNOTATION_COLUMN_ONSET_FLOAT)
  415. eeg_snippet_data[BEHAVIORAL_COLUMN_FIXATION_START] = eeg_snippet_data[EEG_STIMULUS_FIXATION_CROSS]
  416. eeg_snippet_data[BEHAVIORAL_COLUMN_START] = eeg_snippet_data.apply(
  417. lambda row: [row[stimulus] for stimulus in EEG_STIMULUS_SNIPPET_START.values() if not pd.isna(row[stimulus])][0], axis=1)
  418. eeg_snippet_data[BEHAVIORAL_COLUMN_END] = eeg_snippet_data[EEG_STIMULUS_SNIPPET_END]
  419. eeg_snippet_data = eeg_snippet_data[[
  420. BEHAVIORAL_COLUMN_FIXATION_START, BEHAVIORAL_COLUMN_START, BEHAVIORAL_COLUMN_END]]
  421. return eeg_snippet_data
  424. def split_eeg_segments(eeg_data: Raw, trial_annotations: pd.DataFrame) -> dict[str, Raw]:
  425. '''Split given Raw object to receive eeg segments per snippet.
  427. Each split starts shortly before first fixation cross,
  428. and end shortly after next snippet end, or until the next fixation cross.
  429. If there is no other way of determination, a long buffer is added to the current fixation crops of the split.
  431. Arguments:
  432. * eeg_data: the eeg data
  433. * sequence_order: the snippet sequence order to assign eeg segments to each snippet
  434. * time_before_first_snippet: the time buffer to add before the first fixation cross
  435. * time_after_last_snippet: the time buffer to add after the last snippet end if found
  436. * time_constant_without_ending: the time buffer to add after the last fixation cross if no end found
  438. returns: eeg split per snippet
  439. '''
  440. eeg_segments = {}
  441. for snippet, row in trial_annotations.iterrows():
  442. start = 0
  443. if not pd.isna(row[BEHAVIORAL_COLUMN_FIXATION_START]):
  444. start = row[BEHAVIORAL_COLUMN_FIXATION_START]
  445. elif not pd.isna(row[BEHAVIORAL_COLUMN_START]):
  446. start = row[BEHAVIORAL_COLUMN_START]-5.0
  447. else:
  448. print(
  449. f'ignored snippet {snippet} of row {row} due to missing start')
  450. continue
  451. start = max(0, start)
  452. end = 0
  453. if not pd.isna(row[BEHAVIORAL_COLUMN_END]):
  454. end = row[BEHAVIORAL_COLUMN_END]
  455. else:
  456. end = start+EEG_LONG_BUFFER
  457. end = min(end, eeg_data.times[-1])
  458. if start >= end:
  459. print(start, end)
  460. eeg_segment_data: Raw = eeg_data.copy().crop(
  461. tmin=start, tmax=end, include_tmax=True)
  462. # remove unneeded annotations
  463. segment_annotations = eeg_segment_data.annotations.to_data_frame()
  464. segment_annotations['unnecessary'] = ~segment_annotations[ANNOTATION_COLUMN_DESCRIPTION].isin(
  465. list(EEG_STIMULUS_SNIPPET_START.values())+[EEG_STIMULUS_FIXATION_CROSS, EEG_STIMULUS_SNIPPET_END])
  466. eeg_segment_data.annotations.delete(
  467. segment_annotations[segment_annotations['unnecessary']].index)
  468. # add begin annotation to add snippet name to file
  469. eeg_segment_data.annotations.append(
  470. start, 1/EEG_FREQUENCY, f'SNIPPET {snippet}')
  471. eeg_segments[snippet] = eeg_segment_data
  472. return eeg_segments
  475. def check_voltage_per_segment(eeg_segments: dict[str, Raw]) -> pd.DataFrame:
  476. snippet_violation_data = pd.DataFrame(index=[snippet for snippet in eeg_segments],
  477. columns=['Voltage Step Count', 'Voltage Step Channels',
  478. 'Voltage Step Frames', 'Voltage Difference Count',
  479. 'Voltage Difference Channels', 'Voltage Difference Frames'], dtype=object)
  481. for snippet in eeg_segments:
  482. eeg_segment = eeg_segments[snippet]
  483. snippet_violations = check_voltage_in_segment(eeg_segment)
  484. for key in snippet_violations:
  485. snippet_violation_data.at[snippet, key] = snippet_violations[key]
  486. return snippet_violation_data
  489. def check_voltage_in_segment(eeg_segment: Raw, is_epoch: bool = False) -> pd.DataFrame:
  490. snippet_violations = {}
  491. assert eeg_segment.info[MNE_KEY_FREQUENCY] == EEG_FREQUENCY
  492. # check voltage steps
  493. eeg_content_data: np.ndarray = eeg_segment.get_data(
  494. picks=EEG_CHANNELS,units='uV')
  495. if is_epoch:
  496. eeg_content_data = eeg_content_data[0]
  497. # * voltage steps >= 30µV/1ms --> (or 60µV/2ms ?)
  498. voltage_step = np.diff(eeg_content_data, 1)
  499. abs_voltage_step = np.abs(voltage_step)
  500. high_voltage_step = abs_voltage_step >= EEG_VOLTAGE_STEP
  501. has_high_voltage_step = np.any(high_voltage_step)
  502. if has_high_voltage_step:
  503. count_high_voltage_step = np.sum(high_voltage_step*1)
  504. snippet_violations['Voltage Step Count'] = count_high_voltage_step
  505. channel_high_voltage_step = np.sum(high_voltage_step*1, -1)
  506. snippet_violations['Voltage Step Channels'] = [
  507. {channel: channel_high_voltage_step[i] for i, channel in enumerate(EEG_CHANNELS) if channel_high_voltage_step[i] > 0}]
  508. time_high_voltage_step = np.sum(high_voltage_step*1, 0)
  509. frame_high_voltage_step = np.nonzero(time_high_voltage_step)[0]
  510. snippet_violations['Voltage Step Frames'] = [
  511. {frame: time_high_voltage_step[frame] for frame in frame_high_voltage_step}]
  512. starts = [eeg_segment.first_time + f /
  513. EEG_FREQUENCY for f in frame_high_voltage_step]
  514. eeg_segment.annotations.append(
  515. starts, 1/EEG_FREQUENCY, 'BAD Voltage step')
  516. # plot_eeg(eeg_segment, True, True)
  517. del count_high_voltage_step
  518. del time_high_voltage_step
  519. if is_epoch:
  520. return snippet_violations
  521. del voltage_step
  522. del abs_voltage_step
  523. del high_voltage_step
  524. # * voltage difference > 100 µV within 0.2 s
  525. voltage_windows = np.lib.stride_tricks.sliding_window_view(
  526. eeg_content_data, 101, -1)
  527. min_voltage_windows = np.min(voltage_windows, 2)
  528. max_voltage_windows = np.max(voltage_windows, 2)
  529. high_difference_voltage = (
  530. max_voltage_windows-min_voltage_windows) > EEG_VOLTAGE_WINDOW
  531. has_high_difference_voltage = np.any(high_difference_voltage)
  532. if has_high_difference_voltage:
  533. count_high_difference_voltage = np.sum(high_difference_voltage*1)
  534. snippet_violations['Voltage Difference Count'] = count_high_difference_voltage
  535. channel_high_difference_voltage = np.sum(
  536. high_difference_voltage*1, -1)
  537. snippet_violations['Voltage Difference Channels'] = [
  538. {channel: channel_high_difference_voltage[i] for i, channel in enumerate(EEG_CHANNELS) if channel_high_difference_voltage[i] > 0}]
  539. time_high_difference_voltage = np.sum(high_difference_voltage*1, 0)
  540. frame_high_difference_voltage = np.nonzero(
  541. time_high_difference_voltage)[0]
  542. # snippet_violation_data.loc[snippet,'Voltage Difference Frames'] =[{frame:time_high_difference_voltage[frame] for frame in frame_high_difference_voltage.flat}]
  543. first_time = eeg_segment.tmin if is_epoch else eeg_segment.first_time
  544. starts = [first_time + f /
  545. EEG_FREQUENCY for f in frame_high_difference_voltage]
  546. last_time = eeg_segment.tmax if is_epoch else eeg_segment._last_time
  547. durations = [
  548. min(s+100/EEG_FREQUENCY, last_time)-s for s in starts]
  549. eeg_segment.annotations.append(
  550. starts, durations, 'BAD Voltage difference')
  551. # print(eeg_segment.annotations.to_data_frame())
  552. # plot_eeg(eeg_segment, True, True)
  553. del count_high_difference_voltage
  554. del time_high_difference_voltage
  555. del voltage_windows
  556. del min_voltage_windows
  557. del max_voltage_windows
  558. del high_difference_voltage
  559. return snippet_violations
  562. def check_voltage_amplitude(epochs) -> bool:
  563. # * greater absolute amplitude difference than 140 µV --> or if baseline corrected, then within +/-70 µV
  564. eeg_content_data = epochs[0].get_data(picks=EEG_CHANNELS,units='uV')[0]
  565. min_overall_voltage = np.min(eeg_content_data, (0, 1))
  566. max_overall_voltage = np.max(eeg_content_data, (0, 1))
  567. if (min_overall_voltage < EEG_VOLTAGE_OVERALL[0]) or (max_overall_voltage > EEG_VOLTAGE_OVERALL[1]):
  568. return True
  569. return False
  572. def perform_eeg_erp_averaging(participants: list[str], erp_frp: bool | str = True, correct_data_only: bool = False, epoch_interval: tuple[int, int] = (-0.2, 1),
  573. conditional_stimuli: dict[str, str] = EEG_STIMULUS_SNIPPET_START, topomap_times: list[float] = [0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1],
  574. plot: bool = False, snippet_group: str = SNIPPET_GROUP_ALL, snippet_numbers: list[int] = SNIPPET_NUMBERS) -> tuple[dict[str, dict[str, mne.Evoked], dict[str, mne.Evoked]]]:
  575. assert (erp_frp is True or erp_frp in FIXATION_SELECTION_ALGORITHMS), erp_frp
  576. assert epoch_interval[0] < epoch_interval[1]
  577. assert all([t >= epoch_interval[0] and t <= epoch_interval[1]
  578. for t in topomap_times])
  579. description = get_erp_description(
  580. erp_frp, correct_data_only, epoch_interval)
  581. subjectwise_nave = {}
  582. subjectwise_averages = {}
  583. subjectwise_frp_offsets = {}
  584. for participant in tqdm(participants):
  585. print('----------------------------------------------')
  586. print(PARTICIPANT, participant)
  587. # skip participant if excluded
  588. exclusions = get_exclusions(participant, [PARTICIPANT], [
  589. BEHAVIORAL, EEG, VISUAL])[PARTICIPANT]
  590. if any(exclusions.values()):
  591. print('excluded')
  592. continue
  593. # Load behavioral data
  594. if correct_data_only:
  595. behavioral_data = pd.read_csv(get_behavioral_data_path(
  596. participant, final_data_exclusion=True), index_col=False, sep=SEPARATOR, dtype={PARTICIPANT: str})
  597. correct_snippets = behavioral_data[behavioral_data[BEHAVIORAL_COLUMN_CORRECTNESS]][SNIPPET].to_list()
  598. else:
  599. correct_snippets = None
  600. # Load raw data for all snippets
  601. snippet_segments = get_all_eeg_trial_segments(
  602. erp_frp == True, participant, correct_trials=correct_snippets, accepted_snippet_numbers=snippet_numbers)
  604. # Transform annotations to event to epoch and apply baseline correction
  605. snippet_epochs, frp_fixation_offsets = extract_epochs_from_snippet_segments(erp_frp, description, participant, snippet_group,
  606. snippet_segments, conditional_stimuli.values(), epoch_interval, (epoch_interval[0], 0), True, True)
  607. frp_fixation_offsets[PARTICIPANT] = participant
  608. subjectwise_frp_offsets[participant] = frp_fixation_offsets
  609. # Concatenate epochs
  610. snippet_groups = [list(snippet_epochs.keys())]
  611. group_epochs = mne.concatenate_epochs([snippet_epochs[snippet] for snippet in snippet_groups[0]], add_offset=True).pick(
  612. picks=EEG_CHANNELS).set_montage("easycap-M1")
  613. # Calculate and plot subjectwise average per condition
  614. averaged_evoked: mne.Evoked = group_epochs.average(by_event_type=True)
  615. averaged_evoked = {ev.comment: ev for ev in averaged_evoked}
  616. averaged_evoked = {
  617. condition: averaged_evoked[conditional_stimuli[condition]] for condition in conditional_stimuli}
  618. assert (averaged_evoked[CONDITION_CLEAN].comment == conditional_stimuli[CONDITION_CLEAN]) and (averaged_evoked[CONDITION_CONFUSING].comment == conditional_stimuli[CONDITION_CONFUSING]), \
  619. f'''{averaged_evoked[CONDITION_CLEAN].comment} should be {conditional_stimuli[CONDITION_CLEAN]} and {
  620. averaged_evoked[CONDITION_CONFUSING].comment} should be {conditional_stimuli[CONDITION_CONFUSING]}'''
  621. subjectwise_averages[participant] = averaged_evoked
  622. subjectwise_nave[participant] = {
  623. condition: averaged_evoked[condition].nave for condition in averaged_evoked}
  624. if plot:
  625. plot_all_evoked_low_frequency(erp_frp, description, averaged_evoked,
  626. participant, topomap_times, False)
  627. # comparison to BrainVision results
  628. # ae_c,_ = load_eeg_data(f'E:/PHD/Studies/aoc-frp-main-studies/Main_Study_Part/08-Data-Trial_Recordings/prepared_EEG_files/ERP subjectwise averages/AoCfrp_{participant}__averaged_confusing.vhdr', preload=True)
  629. # ae_c_events, ae_c_event_dict = mne.events_from_annotations(ae_c, {'Time 0/':11})
  630. # ae_c_epochs = mne.Epochs(ae_c, ae_c_events, tmin=ERP_INTERVAL[0]+1/EEG_FREQUENCY, tmax=ERP_INTERVAL[1]-1/EEG_FREQUENCY, event_id=ae_c_event_dict, preload=True, baseline=(ERP_INTERVAL[0]+1/EEG_FREQUENCY,0))
  631. # ae_c_averaged_epochs =ae_c_epochs.average()
  632. # plot_evoked(ae_c_averaged_epochs)
  633. # ae_nc,_ = load_eeg_data(f'E:/PHD/Studies/aoc-frp-main-studies/Main_Study_Part/08-Data-Trial_Recordings/prepared_EEG_files/ERP subjectwise averages/AoCfrp_{participant}__averaged_non_confusing.vhdr', preload=True)
  634. # ae_nc_events, ae_nc_event_dict = mne.events_from_annotations(ae_nc, {'Time 0/':12})
  635. # ae_nc_epochs = mne.Epochs(ae_nc, ae_nc_events, tmin=ERP_INTERVAL[0]+1/EEG_FREQUENCY, tmax=ERP_INTERVAL[1]-1/EEG_FREQUENCY, event_id=ae_nc_event_dict, preload=True, baseline=(ERP_INTERVAL[0]+1/EEG_FREQUENCY,0))
  636. # ae_nc_averaged_epochs =ae_nc_epochs.average()
  637. # plot_evoked(ae_nc_averaged_epochs)
  638. # Calculate and plot subjectwise difference wave
  639. diff_wave = mne.combine_evoked(
  640. [averaged_evoked[CONDITION_CONFUSING], averaged_evoked[CONDITION_CLEAN]], weights=[1, -1])
  641. subjectwise_averages[participant][CONDITION_DIFF] = diff_wave
  642. if plot:
  643. plot_all_evoked_low_frequency(erp_frp, description, {
  644. CONDITION_DIFF: diff_wave}, participant, topomap_times, False)
  645. # Save subjectwise averages
  646. for condition in averaged_evoked:
  647. averaged_evoked[condition].save(get_erp_average_path(
  648. erp_frp, snippet_group, description, participant, condition), overwrite=True)
  649. # save included snippets, offsets
  650. frp_offset_data = pd.concat(subjectwise_frp_offsets.values())
  651. frp_offset_data[ERP_PARAMETER_EPOCH_INTERVAL] = [
  652. epoch_interval for _ in range(frp_offset_data.shape[0])]
  653. frp_offset_data[ERP_PARAMETER_CORRECT_TRIALS_ONLY] = correct_data_only
  654. if erp_frp != True:
  655. frp_offset_data[FIXATION_SELECTION_ALGORITHM] = erp_frp
  656. frp_offset_data.to_csv(get_erp_fixation_analysis_path(
  657. erp_frp, snippet_group, description, 'erp frp offset'), sep=SEPARATOR, index=False)
  658. # statistics and plot distribution (best in other method)
  659. if erp_frp != True:
  660. statistics_distribution(erp_frp, snippet_group, description,
  661. frp_offset_data, 'erp frp offset', 'Delay to stimulus onset')
  662. # Save subjectwise naves
  663. nave_data = pd.DataFrame.from_dict(subjectwise_nave, 'index')
  664. nave_data.to_csv(get_erp_nave_path(
  665. erp_frp, snippet_group, description), sep=SEPARATOR)
  666. # Calculate and plot grand averages per condition
  667. grand_averages = {}
  668. for condition in [CONDITION_CONFUSING, CONDITION_CLEAN]:
  669. grand_average = mne.grand_average(
  670. [subjectwise_averages[participant][condition] for participant in subjectwise_averages])
  671. grand_averages[condition] = grand_average
  672. if plot:
  673. plot_all_evoked_low_frequency(erp_frp, description, grand_averages,
  674. TOTAL, topomap_times, False)
  675. # Calculate and plot grand averages difference wave
  676. diff_wave = mne.combine_evoked(
  677. [grand_averages[CONDITION_CONFUSING], grand_averages[CONDITION_CLEAN]], weights=[1, -1])
  678. grand_averages[CONDITION_DIFF] = diff_wave
  679. if plot:
  680. plot_all_evoked_low_frequency(erp_frp, description, {
  681. CONDITION_DIFF: diff_wave}, TOTAL, topomap_times, False)
  682. # Save grand averages
  683. for condition in grand_averages:
  684. grand_averages[condition].save(get_erp_average_path(
  685. erp_frp, snippet_group, description, TOTAL, condition=condition), overwrite=True)
  686. return subjectwise_averages, grand_averages, subjectwise_nave
  689. def get_erp_description(erp_frp: str, correct_data_only: bool, epoch_interval: tuple[int, int]) -> str:
  690. description = f'{"erp" if erp_frp is True else FIXATION_SELECTION_SHORT_VERSION[erp_frp]}_{int(epoch_interval[0]*1000)}_{int(epoch_interval[1]*1000)}_{"correct" if correct_data_only else "all"}'
  691. return description
  694. def get_stimulus_number(stimuli: list[str] = EEG_STIMULUS) -> dict[str, int]:
  695. '''get stimulus number (event number) per recognized stimulus (used in description of annotations)
  697. returns: stimulus text and number per recognized stimulus'''
  698. return {stimulus: int(stimulus[-3:]) for stimulus in stimuli}
  701. def get_event_name_numbers(given_stimuli: list[str]) -> dict[str, int]:
  702. '''get event name per event number for each recognized stimulus (used in description of annotations)
  704. Arguments: given_stimuli: the stimuli to return
  706. returns: event name and event number per recognized stimulus'''
  707. given_stimuli = [
  708. stimulus for stimulus in given_stimuli if stimulus in EEG_STIMULUS]
  709. stimulus_numbers = get_stimulus_number()
  710. return {STIMULUS_EVENT_NAMES[stimulus]: stimulus_numbers[stimulus] for stimulus in given_stimuli}
  713. def get_all_eeg_trial_segments(erp_frp: bool | str, participant: str, correct_trials: list[str] = None, accepted_snippet_numbers: list[int] = None) -> dict[str, Raw]:
  714. '''load eeg trial segments for participant
  716. Arguments:
  717. * erp_frp: whether it is erp or a certain type of frp
  718. * participant: the participant to get segments for
  719. * visual_exclude: whether to exclude based on visual as well (required for FRP)
  721. returns: the non-excluded trial segments for this participant
  722. '''
  723. assert (erp_frp in [True, False]), erp_frp
  724. modes = [BEHAVIORAL, EEG]
  725. if not erp_frp:
  726. modes.append(VISUAL)
  727. exclusions = get_exclusions(participant, [SNIPPET], modes)[SNIPPET]
  728. snippet_segments: dict[str, Raw] = {}
  729. for snippet in exclusions:
  730. if any(exclusions[snippet].values()):
  731. # print(f'\t{snippet} was excluded due to {exclusions[snippet]}')
  732. continue
  733. # if correct only (correct trials given) and snippet not correctly answered
  734. if not (correct_trials is None) and not (snippet in correct_trials):
  735. # print(f'\t{snippet} was excluded due to being answered incorrectly')
  736. continue
  737. # if correct only (correct trials given) and snippet not correctly answered
  738. if not (accepted_snippet_numbers is None) and not (get_snippet_number(snippet) in accepted_snippet_numbers):
  739. # print(f'\t{snippet} was excluded due to not being in the group')
  740. continue
  741. eeg_data, _ = load_eeg_data(get_eeg_trial_path(
  742. erp_frp is True, participant, snippet))
  743. snippet_segments[snippet] = eeg_data
  744. return snippet_segments
  747. def plot_epoch(eeg_data: mne.Epochs):
  748. return eeg_data.plot(EEG_CHANNELS, n_epochs=1, events=eeg_data.events)
  751. def round_time_EEG(time: float, frequency=EEG_FREQUENCY):
  752. return round(time, int(round(log10(frequency), 0))+1)
  755. def extract_epochs_from_snippet_segments(erp_frp: bool | str, description: str, participant: str, snippet_group: str, snippet_segments: dict[str, Raw], condition_stimuli: list[str],
  756. epoch_interval: tuple[int, int], baseline_interval: tuple[int, int] = None, perform_voltage_checks: bool = True, save_epoch_data: bool = True):
  757. regarded_stimuli = get_stimulus_number(condition_stimuli)
  758. snippet_epochs = {}
  759. frp_fixation_offsets = pd.DataFrame([], columns=[
  760. PARTICIPANT, CONDITION, SNIPPET, 'Stimulus Onset', 'Fixation Onset', 'Delay to stimulus onset'])
  761. snippet_status = {}
  762. # calculate epoch
  763. for snippet in snippet_segments:
  764. try:
  765. events, event_dict = mne.events_from_annotations(
  766. snippet_segments[snippet], regarded_stimuli)
  767. except ValueError:
  768. if erp_frp != True:
  769. print(
  770. f'''{snippet}: No stimulus of {list(condition_stimuli)} found in annotations in {snippet_segments[snippet].annotations.to_data_frame()[ANNOTATION_COLUMN_DESCRIPTION].values}''')
  771. snippet_status[snippet] = 'No stimulus found, fixation data of this trials likely did not contain any fixation fulfilling the requirements for this FRP calculation.'
  772. continue
  773. else:
  774. raise ValueError(
  775. f'''In ERP, all stimuli must be found. No stimulus of {list(condition_stimuli)} found in annotations for {snippet} {snippet_segments[snippet].annotations.to_data_frame()[ANNOTATION_COLUMN_DESCRIPTION].values}''')
  776. # create epoch and perform baseline correction if specified
  777. if baseline_interval is None:
  778. epochs = mne.Epochs(snippet_segments[snippet], events, tmin=epoch_interval[0],
  779. tmax=epoch_interval[1], event_id=event_dict, preload=True)
  780. else:
  781. epochs = mne.Epochs(snippet_segments[snippet], events, tmin=epoch_interval[0], tmax=epoch_interval[1],
  782. event_id=event_dict, preload=True, baseline=(baseline_interval[0], baseline_interval[1]))
  783. # check whether epoch really exists (not too short)
  784. if erp_frp != True and len(epochs) < 1 and 'TOO_SHORT' in epochs.drop_log[0]:
  785. anno = snippet_segments[snippet].annotations.to_data_frame()
  786. relevant_stimuli = anno[anno[ANNOTATION_COLUMN_DESCRIPTION].isin([EEG_STIMULUS_SNIPPET_END, *EEG_STIMULUS_SNIPPET_START.values(), *condition_stimuli])]
  787. snippet_status[snippet] = f'Data for existing stimuli too short, {relevant_stimuli}'
  788. print(
  789. f'''{snippet}: Data for existing stimuli {list(event_dict.keys())[0]} too short {epochs.drop_log[0]} {relevant_stimuli}''')
  790. continue
  791. assert len(epochs) == 1
  792. # exclude based on all previously marked violations
  793. if perform_voltage_checks:
  794. if check_voltage_amplitude(epochs):
  795. print(f'\t{snippet} excluded due to overall voltage violation')
  796. # plot_epoch(epochs)
  797. snippet_segments[snippet].close()
  798. snippet_status[snippet] = 'Voltage violation absolute of segment'
  799. continue
  800. if check_voltage_in_segment(epochs, True):
  801. print(
  802. f'\t{snippet} excluded due to voltage violation inside epoch')
  803. # plot_epoch(epochs)+-
  804. snippet_segments[snippet].close()
  805. snippet_status[snippet] = 'Voltage violation in interval of segment'
  806. continue
  807. snippet_status[snippet] = 'Included'
  808. snippet_epochs[snippet] = epochs
  809. # calculate frp offset to erp
  810. annotation_data, _ = prepare_annotation_information(
  811. snippet_segments[snippet])
  812. erp_onset = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin(
  813. EEG_STIMULUS_SNIPPET_START.values())][ANNOTATION_COLUMN_ONSET_FLOAT].values[0]
  814. if erp_frp != True:
  815. frp_onset = annotation_data[annotation_data[ANNOTATION_COLUMN_DESCRIPTION].isin(
  816. condition_stimuli)][ANNOTATION_COLUMN_ONSET_FLOAT].values[0]
  817. frp_fixation_offsets.loc[frp_fixation_offsets.shape[0]] = [None, CONDITION_VARIANT_MATCH[get_snippet_variant(
  818. snippet)], snippet, round_time_EEG(erp_onset), round_time_EEG(frp_onset), round_time_EEG(frp_onset-erp_onset)]
  819. else:
  820. frp_fixation_offsets.loc[frp_fixation_offsets.shape[0]] = [None, CONDITION_VARIANT_MATCH[get_snippet_variant(
  821. snippet)], snippet, round_time_EEG(erp_onset), round_time_EEG(erp_onset), .0]
  822. with open(get_erp_status_path(erp_frp, snippet_group, description, participant), 'w') as f:
  823. json.dump(snippet_status, f, indent=4, sort_keys=True)
  824. if save_epoch_data:
  825. # save all epochs
  826. for snippet, epoch in snippet_epochs.items():
  827. epoch.save(get_erp_epoch_path(erp_frp, snippet_group, description,
  828. participant, snippet), fmt='double', overwrite=True)
  829. return snippet_epochs, frp_fixation_offsets
  832. def plot_all_evoked(erp_frp: bool | str, snippet_group, description: str, conditional_evoked: dict[str, mne.Evoked], participant: str = TOTAL, topomap_times: list[float] = [0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1], show: bool = True):
  833. assert (erp_frp is True or erp_frp in FIXATION_SELECTION_ALGORITHMS), erp_frp
  834. for condition, evoked in conditional_evoked.items():
  835. fig = evoked.plot(picks='eeg', show=show,
  836. window_title=condition, time_unit='ms')
  837. fig.savefig(get_erp_average_path(erp_frp, snippet_group, description,
  838. participant, condition, 'butterfly'))
  839. plt.close()
  840. for condition, evoked in conditional_evoked.items():
  841. fig1 = evoked.plot_topomap(times=[min(
  842. time, evoked.tmax) for time in topomap_times], time_unit='ms', show=False)
  843. # TODO: hier Daten für Topoplots abgreifen
  844. fig1.suptitle(f'Topomap {description}')
  845. fig1.savefig(get_erp_average_path(erp_frp, snippet_group, description,
  846. participant, condition, 'topomap'))
  847. fig1.show()
  848. if erp_frp!=True:
  849. fig2 = evoked.plot_topomap(times=[min(
  850. round_time_EEG(time+0.0255, 10000), evoked.tmax) for time in topomap_times],average=0.05199, time_unit='ms', show=False)
  851. fig2_path = get_erp_average_path(erp_frp, snippet_group, description,
  852. participant, condition, 'topomap_averaged_50ms')
  853. for i, ax in enumerate(fig2.get_axes()[:-1]):
  854. title:str = ax.get_title()[:-3]
  855. start, end = [int(i)/1000 for i in title.split(' – ')]
  856. new_title = fig1.axes[i].get_title()
  857. data = evoked.copy().crop(tmin=start, tmax=end).to_data_frame(index='time', time_format='ms')
  858. data.index.name='Time (ms)'
  859. data.to_csv(fig2_path.with_name(f'Data Figure2b amplitudes {new_title} interval.csv'), sep=SEPARATOR)
  860. ax.set_title(new_title)
  861. fig2.suptitle(f'Topomap {description}')
  862. fig2.savefig(fig2_path)
  863. fig2.show()
  864. else:
  865. fig2 = evoked.plot_topomap(times=[min(
  866. round_time_EEG(time+0.1015, 10000), evoked.tmax) for time in topomap_times],average=0.201, time_unit='ms', show=False)
  867. fig2_path = get_erp_average_path(erp_frp, snippet_group, description,
  868. participant, condition, 'topomap_averaged_200ms')
  869. for i, ax in enumerate(fig2.get_axes()[:-1]):
  870. title:str = ax.get_title()[:-3]
  871. start, end = [int(i)/1000 for i in title.split(' – ')]
  872. new_title = fig1.axes[i].get_title()
  873. data = evoked.copy().crop(tmin=start, tmax=end).to_data_frame(index='time', time_format='ms')
  874. data.index.name='Time (ms)'
  875. data.to_csv(fig2_path.with_name(f'Data Figure3b amplitudes {new_title} interval.csv'), sep=SEPARATOR)
  876. ax.set_title(new_title)
  877. fig2.suptitle(f'Topomap {description}')
  878. fig2.savefig(fig2_path)
  879. fig2.show()
  880. plt.close('all')
  881. if all([condition in CONDITION_COLORS for condition in conditional_evoked]):
  882. colors = {condition: CONDITION_COLORS[condition]
  883. for condition in conditional_evoked}
  884. else:
  885. colors = None
  886. fig = mne.viz.plot_compare_evokeds(
  887. conditional_evoked, show_sensors=True, title=f'Topographic comparison {description}', axes='topo', show=show, colors=colors, time_unit='ms')
  888. fig[0].savefig(get_erp_average_path(erp_frp, snippet_group, description, participant, (condition if len(
  889. conditional_evoked) == 1 else TOTAL), f'topo_channels'))
  890. plt.close()
  891. minimum, maximum = [], []
  892. for condition, evoked in conditional_evoked.items():
  893. eeg_data = evoked.get_data(EEG_CHANNELS,units='uV')
  894. minimum.append(np.min(eeg_data))
  895. maximum.append(np.max(eeg_data))
  896. minimum, maximum = min(minimum), max(maximum) # from volt to microvolt scale
  897. for channel in EEG_CHANNELS:
  898. fig = mne.viz.plot_compare_evokeds(conditional_evoked, picks=channel, title=f'Electrode {channel}', show=show, colors=colors,
  899. show_sensors=False,
  900. # for identical scaling
  901. ylim={'eeg': (minimum, maximum)}, time_unit='ms', )
  902. fig[0].get_axes()[0].get_legend().remove()
  903. fig_path = get_erp_average_path(erp_frp, snippet_group, description, participant, (condition if len(
  904. conditional_evoked) == 1 else TOTAL), f'channel_{channel}')
  905. fig[0].savefig(fig_path)
  906. plt.close()
  907. if CONDITION_DIFF in conditional_evoked:
  908. return
  909. fig_data_path = get_erp_average_path(erp_frp, snippet_group, description, participant, (condition if len(
  910. conditional_evoked) == 1 else TOTAL), f'data').with_suffix('.csv')
  911. for channel in EEG_CHANNELS:
  912. if channel[0] in 'FCP' and channel[1] in '34z':
  913. data = []
  914. for condition in conditional_evoked:
  915. cond_data = conditional_evoked[condition].to_data_frame(channel, index='time', time_format='ms')
  916. cond_data.columns = [condition]
  917. data.append(cond_data)
  918. data = pd.concat(data, axis=1)
  919. data.index.name = 'Time (ms)'
  920. data.to_csv(fig_data_path.with_stem(f'Data Figure{3 if erp_frp==True else 2}a {channel} conditional amplitude'), sep=SEPARATOR)
  922. def plot_all_evoked_low_frequency(erp_frp: bool | str, snippet_group: str, description: str, conditional_evoked: dict[str, mne.Evoked], participant: str = TOTAL, topomap_times: list[float] = [0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1], show: bool = False):
  923. # plot_all_evoked(erp_frp, description, conditional_evoked, participant, topomap_times, show)
  924. conditional_evoked = {condition: evoked.copy().resample(
  925. 20) for condition, evoked in conditional_evoked.items()}
  926. plot_all_evoked(erp_frp, snippet_group,
  927. f'{description}_20Hz', conditional_evoked, participant, topomap_times, show)
  930. def statistics_distribution(erp_frp: bool | str, snippet_group: str, description: str, fixation_analysis_data: pd.DataFrame, analysis_topic: str, analysis_column: str):
  931. fixation_analysis_data[f'{analysis_column} (ms)'] = fixation_analysis_data[analysis_column]*1000
  932. analysis_column = f'{analysis_column} (ms)'
  933. conditional_offset = fixation_analysis_data[[CONDITION, analysis_column]].groupby(
  934. [CONDITION]).agg({analysis_column: PANDAS_DESCRIPTION_AGG_FUNCTIONS})
  935. conditional_offset.columns = PANDAS_DESCRIPTION_AGG_NAMES
  936. conditional_offset.to_csv(get_erp_fixation_analysis_path(
  937. erp_frp, snippet_group, f'{description}_statistics', analysis_topic), sep=SEPARATOR, decimal=',')
  938. fig, axis = plt.subplots(1, 1, figsize=(8, 3))
  939. plt.rcParams.update({'font.size': 12})
  940. # , palette=[CONDITION_COLORS[CONDITION_CLEAN], CONDITION_COLORS[CONDITION_CONFUSING]])
  941. sns.violinplot(fixation_analysis_data, x=analysis_column, y=CONDITION, legend=False, inner="box", cut=0, ax=axis)
  942. plt.tight_layout()
  943. plt.savefig(get_erp_fixation_analysis_path(erp_frp, snippet_group,
  944. f'{description}_statistics', analysis_topic).with_suffix('.pdf'), bbox_inches='tight', pad_inches=0)
  945. plt.savefig(get_erp_fixation_analysis_path(erp_frp, snippet_group,
  946. f'{description}_statistics', analysis_topic).with_suffix('.png'), bbox_inches='tight', pad_inches=0)
  947. plt.close()
  948. participant_conditional_offset = fixation_analysis_data[[PARTICIPANT, CONDITION, analysis_column]].groupby(
  949. [PARTICIPANT, CONDITION]).agg({analysis_column: PANDAS_DESCRIPTION_AGG_FUNCTIONS})
  950. participant_conditional_offset.columns = PANDAS_DESCRIPTION_AGG_NAMES
  951. participant_conditional_offset.to_csv(get_erp_fixation_analysis_path(
  952. erp_frp, snippet_group, f'{description}_participant_statistics', analysis_topic), sep=SEPARATOR, decimal=',')
  953. fig, axis = plt.subplots(1, 1, figsize=(8, 24))
  954. plt.rcParams.update({'font.size': 12})
  955. sns.violinplot(fixation_analysis_data, x=analysis_column,
  956. y=PARTICIPANT, hue=CONDITION, inner="stick", cut=0, ax=axis)
  957. plt.tight_layout()
  958. plt.savefig(get_erp_fixation_analysis_path(erp_frp, snippet_group,
  959. f'{description}_participant_statistics', analysis_topic).with_suffix('.pdf'), bbox_inches='tight', pad_inches=0)
  960. plt.savefig(get_erp_fixation_analysis_path(erp_frp, snippet_group,
  961. f'{description}_participant_statistics', analysis_topic).with_suffix('.png'), bbox_inches='tight', pad_inches=0)
  962. plt.close()
  963. return conditional_offset, participant_conditional_offset
  966. def load_all_erp_averages(erp_frp: bool | str, snippet_group: str, correct_data_only: bool, epoch_interval: tuple[int, int], subjectwise: bool,
  967. participants: list[str], grand: bool, conditional: bool, diff: bool) -> Union[tuple[dict[str, dict[str, mne.Evoked], dict[str, mne.Evoked]]], dict[str, dict[str, mne.Evoked]], dict[str, mne.Evoked]]:
  968. assert (erp_frp is True or erp_frp in FIXATION_SELECTION_ALGORITHMS), erp_frp
  969. assert ((not subjectwise) or (len(participants) > 0)
  970. ), "if subjectwise participants are required, send with the participants to use"
  971. assert (subjectwise or grand), "Subjectwise or grand or both must be chosen"
  972. description = get_erp_description(
  973. erp_frp, correct_data_only, epoch_interval)
  974. conditions = []
  975. if conditional:
  976. conditions.extend([CONDITION_CONFUSING, CONDITION_CLEAN])
  977. if diff:
  978. conditions.append(CONDITION_DIFF)
  979. if grand:
  980. grand_averages = {}
  981. for condition in conditions:
  982. path = get_erp_average_path(
  983. erp_frp, snippet_group, description, condition=condition)
  984. grand_averages[condition] = mne.read_evokeds(path)[0]
  985. if not subjectwise:
  986. return grand_averages
  987. if subjectwise:
  988. subjectwise_averages = {}
  989. for participant in participants:
  990. averages = {}
  991. for condition in conditions:
  992. path = get_erp_average_path(
  993. erp_frp, snippet_group, description, participant, condition=condition)
  994. averages[condition] = mne.read_evokeds(path)[0]
  995. subjectwise_averages[participant] = averages
  996. if not grand:
  997. return subjectwise_averages
  998. return subjectwise_averages, grand_averages
  1001. def load_all_erp_epochs(erp_frp: bool | str, snippet_group: str, description: str, participants: list[str]) -> dict[str, dict[str, mne.Epochs]]:
  1002. assert (erp_frp is True or erp_frp in FIXATION_SELECTION_ALGORITHMS), erp_frp
  1003. subjectwise_epochs = {}
  1004. for participant in participants:
  1005. snippet_epochs = {}
  1006. snippet_epoch_paths = get_all_erp_epoch_paths(
  1007. erp_frp, snippet_group, description, participant)
  1008. for snippet, epoch_path in snippet_epoch_paths.items():
  1009. snippet_epochs[snippet] = mne.read_epochs(
  1010. epoch_path, proj=True, preload=True, verbose=None)
  1011. subjectwise_epochs[participant] = snippet_epochs
  1012. return subjectwise_epochs
  1015. def plot_eeg(eeg_data: Raw, plot_annotations: bool, plot_data: bool) -> None:
  1016. '''plot eeg data or its events.
  1018. Each split starts shortly before first fixation cross,
  1019. and end shortly after next snippet end, or until the next fixation cross.
  1020. If there is no other way of determination, a long buffer is added to the current fixation crops of the split.
  1022. Arguments:
  1023. * eeg_data: the eeg data
  1024. * plot_annotations: whether to plot the annotations as events
  1025. * plot_data: plot the eeg data
  1026. '''
  1027. events, event_id = mne.events_from_annotations(
  1028. eeg_data, get_stimulus_number())
  1029. event_dict = get_event_name_numbers(event_id.keys())
  1030. event_color = {3: 'r', 4: 'b', 11: 'g', 12: 'y', }
  1031. print('Events ID:', event_id, event_dict)
  1032. if plot_data:
  1033. eeg_data.plot(events=events, start=0, duration=30, color='gray', event_color={k: event_color[k] for k in event_color if k in event_id.values()},
  1034. )
  1035. # prepare data
  1036. # for erp_parameters in tqdm(all_erp_parameter_combinations):
  1039. def plot_waveforms(erp_frp: bool | str, snippet_group: str, correct_data_only: bool, epoch_interval: tuple[float, float], subjectwise: bool, grand: bool, participants: list[str], topomap_times: list[float]):
  1040. assert (erp_frp is True or erp_frp in FIXATION_SELECTION_ALGORITHMS), erp_frp
  1041. description = get_erp_description(
  1042. erp_frp, correct_data_only, epoch_interval)
  1043. print(description)
  1044. data = load_all_erp_averages(erp_frp, snippet_group, correct_data_only, epoch_interval,
  1045. subjectwise=subjectwise, participants=participants, grand=grand, diff=True, conditional=True)
  1046. if subjectwise and grand:
  1047. subjectwise_data, total_data = data
  1048. elif subjectwise:
  1049. subjectwise_data = data
  1050. elif grand:
  1051. total_data = data
  1052. if subjectwise:
  1053. for participant in tqdm(participants):
  1054. plot_all_evoked(erp_frp, snippet_group, description, {condition: subjectwise_data[participant][condition] for condition in [CONDITION_CLEAN, CONDITION_CONFUSING]},
  1055. participant, topomap_times, False)
  1056. plot_all_evoked(erp_frp, snippet_group, description, {CONDITION_DIFF: subjectwise_data[participant][CONDITION_DIFF]},
  1057. participant, topomap_times, False)
  1058. del subjectwise_data
  1059. if grand:
  1060. plot_all_evoked(erp_frp, snippet_group, description, {condition: total_data[condition] for condition in [CONDITION_CLEAN, CONDITION_CONFUSING]},
  1061. TOTAL, topomap_times, False)
  1062. plot_all_evoked(erp_frp, snippet_group, description, {CONDITION_DIFF: total_data[CONDITION_DIFF]},
  1063. TOTAL, topomap_times, False)
  1064. for channel in EEG_CHANNELS:
  1065. fig, axis = plt.subplots(1, 1, figsize=(8, 3))
  1066. plt.rcParams.update({'font.size': 12})
  1067. mne.viz.plot_compare_evokeds(total_data, picks=channel, # show_sensors=True,
  1068. title=f'Electrode {channel}', show=False, colors=CONDITION_COLORS, time_unit='ms', axes=axis, truncate_yaxis=False)
  1069. fig.savefig(get_erp_average_path(erp_frp, snippet_group, description, TOTAL, 'all',
  1070. f'channel_{channel}').with_suffix('.pdf'), bbox_inches='tight', pad_inches=0)
  1071. plt.clf()
  1072. plt.close()
  1073. total_data[CONDITION_DIFF].plot_topomap(times=[0.400, 0.450, 0.500, 0.550, 0.600, 0.650, 0.700], show=False, time_unit='ms')
  1074. fig.savefig(get_erp_average_path(erp_frp, snippet_group, description,
  1075. TOTAL, CONDITION_DIFF, 'topomap').with_suffix('.pdf'), bbox_inches='tight', pad_inches=0)
  1076. plt.close()
  1077. total_data.pop(CONDITION_DIFF)
  1078. # plt.figure(figsize=(16, 16))
  1079. fig = mne.viz.plot_compare_evokeds(
  1080. {condition:data.copy().crop(tmax=min(data.tmax, 1.0), include_tmax=True) for condition, data in total_data.items()}, picks=['F3','Fz', 'F4','C3','Cz', 'C4','P3','Pz', 'P4', 'F7','F8'],
  1081. show_sensors=True, title=f'Topographic comparison with 9 crucial electrodes {description}', axes='topo',
  1082. colors={condition:CONDITION_COLORS[condition] for condition in total_data}, time_unit='ms', )
  1083. # fig[0].tight_layout()
  1084. fig[0].savefig(get_erp_average_path(erp_frp, snippet_group, description, TOTAL, 'all', 'topo_9_channels_1sec').with_suffix('.pdf'))#, bbox_inches='tight', pad_inches=0)
  1085. fig = mne.viz.plot_compare_evokeds(
  1086. total_data, picks=['F3','Fz', 'F4','C3','Cz', 'C4','P3','Pz', 'P4', 'F7','F8'],
  1087. show_sensors=True, title=f'Topographic comparison with 9 crucial electrodes {description}', axes='topo', show=False,
  1088. colors={condition:CONDITION_COLORS[condition] for condition in total_data}, time_unit='ms', )
  1089. # fig[0].tight_layout()
  1090. fig[0].savefig(get_erp_average_path(erp_frp, snippet_group, description, TOTAL, 'all', 'topo_9_channels').with_suffix('.pdf'))#, bbox_inches='tight', pad_inches=0)
  1091. plt.close()
  1092. del total_data
  1093. gc.collect()
  1096. def add_frp_marker_by_special_fixations(participant: str, snippets: list[str], behavioral_data: pd.DataFrame, special_fixation_data: pd.DataFrame, eeg_trials: dict[str, mne.io.Raw]) -> float:
  1097. for snippet in snippets:
  1098. snippet_behavioral_start = behavioral_data[behavioral_data[SNIPPET] == snippet].squeeze(
  1099. )[BEHAVIORAL_COLUMN_START]
  1100. eeg_segment_data: mne.io.Raw = eeg_trials[snippet]
  1101. snippet_condition = CONDITION_VARIANT_MATCH[get_snippet_variant(
  1102. snippet)]
  1103. eeg_annotations, _ = prepare_annotation_information(eeg_segment_data)
  1104. eeg_start = eeg_annotations[eeg_annotations[ANNOTATION_COLUMN_DESCRIPTION]
  1105. == EEG_STIMULUS_SNIPPET_START[snippet_condition]].squeeze()[ANNOTATION_COLUMN_ONSET_FLOAT]
  1106. snippet_special_fixation_data: pd.DataFrame = special_fixation_data[
  1107. special_fixation_data[SNIPPET] == snippet]
  1108. for f_a in snippet_special_fixation_data[FIXATION_SELECTION_ALGORITHM].unique():
  1109. fixation_start = snippet_special_fixation_data[snippet_special_fixation_data[FIXATION_SELECTION_ALGORITHM] == f_a].squeeze(
  1110. )[FIXATION_COLUMN_START]
  1111. eeg_fixation_start = transform_eye_to_eeg(fixation_start,
  1112. snippet_behavioral_start, eeg_start)
  1113. eeg_segment_data.annotations.append(eeg_fixation_start, 1/EEG_FREQUENCY,
  1114. FRP_EEG_STIMULUS_SNIPPET_START[f_a][snippet_condition])
  1115. export_eeg_brainvision(
  1116. eeg_segment_data, get_eeg_trial_path(False, participant, snippet))
  1119. # transform eye-tracking timestamp to eeg frame
  1120. def transform_eye_to_eeg(eye_timestamp: float, eye_start_time: float, eeg_start_time: float) -> float:
  1121. '''transforms the eye-tracking timestamp into an eeg frame
  1123. Arguments:
  1124. * eye_timestamp: timestamp of eye-tracking to transform (in seconds with milliseconds floating precision)
  1125. * eye_start_time: timestamp of eye-tracking marking the starting point (in seconds with milliseconds floating precision)
  1126. * eeg_start_time: frame that corresponds to the eye_start timestamp
  1127. * eeg_sampling_rate: the frequency of frames logged in the data (frames per second)
  1129. returns: the frame corresponding to the eye_timestamp
  1130. '''
  1131. eye_offset = eye_timestamp-eye_start_time
  1132. eeg_timestamp = eeg_start_time+eye_offset
  1133. return eeg_timestamp

eeg_helpers.py at commit d84bf83, under CC-BY-4.0 · at the source

Overview

Authors: Annabelle Bergum1, Anna-Maria Maurer1, Norman Peitek1, Regine Bader2, Axel Mecklinger2, Vera Demberg1,3, Janet Siegmund4, Sven Apel1
ORCID iDs: Anna-Maria Maurer
  1. Computer Science, Saarland University, Saarbrücken, Germany
  2. Psychology, Saarland University, Saarbrücken, Germany
  3. Language Science and Technology, Saarland University, Saarbrücken, Germany
  4. Computer Science, University of Technology Chemnitz, Chemnitz, Germany
Institutions: Saarland University (Germany); Chemnitz University of Technology (Germany)
Journal: Scientific reports, volume 16, issue 1, article 16833
Dates: received 23 January 2025; accepted 24 April 2026; published online 1 June 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1038/s41598-026-50946-9 · PMID 42225689 · PMCID PMC13226682 · OpenAlex W4405433354
Open access: gold, a free copy (OpenAlex)
Preprint: arxiv:2412.10099
Status: code verified
Categories: EEG (modality), human (organism), methods / tools (subfield)
Methods: Spectral & time-frequency, Preprocessing, Statistics, Smoothing, state filtering, decompositions, Evoked potentials, Physiology & signal measures
Keywords: Program comprehension, Natural language comprehension, Atoms of confusion, EEG, Fixation-related potentials, Late frontal positivity, Software, Computer science, Psychology, Human behaviour
MeSH: Evoked Potentials*, Software*, Comprehension, Electroencephalography, Humans, Natural Language Processing (* major topic)
Topic: Software Engineering Research (Information Systems, Computer Science), according to OpenAlex
Funding: European Research Council (101052182)
Citations: not cited yet (Europe PMC); 70 references in the paper

Abstract

As software pervades more and more areas of our professional and personal lives, there is an ever-increasing need to maintain software and for programmers to efficiently write and understand program code. In the first study of its kind, we analyze fixation-related potentials (FRPs) to explore the online processing of program code patterns that are confusing to programmers, but not to the computer (so-called atoms of confusion), and their underlying neurocognitive mechanisms in an ecologically valid setting. Relative to clean counterparts in program code without an atom of confusion, confusing code elicits a late frontal positivity of about 400 to 700 ms after first looking at the atom of confusion. This frontal positivity resembles an event-related potential (ERP) component found during natural language processing that is elicited by unexpected but plausible words in sentence context. Thus, we suggest that the brain engages similar neurocognitive mechanisms in response to unexpected and informative inputs in program code and in natural language. In both domains, these inputs update a comprehender’s situation model, which is essential for information extraction from a quickly unfolding input. Our results have far-reaching implications for programming and pave the way for interdisciplinary collaborations between software engineering and psycholinguistics.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repositories

Its files are read in the Code ↔ Paper reader above, with 5 matches between paragraphs and lines of code.

Zenodo 14229848

License: CC-BY-4.0
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Size: 8 files
Software Heritage: not checked
Found in: “Data availability”
Holds: README
Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation
Availability: 1 check, the latest on 27 September 2026: the link answers (HTTP 200)
  • 27 September 2026: the link answers (HTTP 200)
1 file
At the source:

brains-on-code/AoC-FRP-Code

License: CC-BY-4.0
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: d84bf83e386dc631437efe80e7d5ac75098a51c8, 8 January 2025
Languages: Java (145), Jupyter (39), Python (26)
Size: 611 files, 210 scripts
Software Heritage: not checked
Found in: “Materials availability”
Holds: README, license file, CITATION.cff, environment (requirements.txt), 39 notebooks
Not found: tests, continuous integration, documentation
Tools: pandas (30 files), NumPy (17 files), Matplotlib (9 files), SciPy (9 files), MNE-Python (7 files), seaborn (3 files), h5py (2 files), OpenCV (2 files), lme4 (1 file), lmerTest (1 file), nlme (1 file), Pillow (1 file), PsychoPy (1 file), rpy2 (1 file), scikit-learn (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
212 files
At the source: github.com/brains-on-code/AoC-FRP-Code

Materials availability

This experiment used code available in Python (v. 3.11.5), as well as R (v. 4.3.2) for linear mixed effect models. Additionally, we used other open-source (PsychoPy, v. 2021.2.3) and commercial applications (Tobii Eye-Tracker Manager (v. 2.6.0), BrainVision Recorder (v. 1.20.0801), and BrainVision Analyzer (v. 2.3.0.8300)) to perform the experiment. We deposited the scripts and content in a GitHub Repository https://github.com/brains-on-code/AoC-FRP-Code and they are available for unrestricted open access under a CC-BY license.

Reproduced under the paper's license (CC BY), from the paper cited above.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 2 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 210 scripts, each with its path and the digest of its content;
  • 5 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

Datasets cited

Data Availability Statement

Data availability We provide all relevant data in line with open data principles under a CC-BY license, respecting our local privacy laws (GDPR). Specifically, the datasets collected during the experiment and generated for the analysis are long-term archived in the Zenodo repository Dataset for “Fixation-related potentials reveal that confusing program code elicits a late frontal positivity”: 10.5281/zenodo.14229848

This experiment used code available in Python (v. 3.11.5), as well as R (v. 4.3.2) for linear mixed effect models. Additionally, we used other open-source (PsychoPy, v. 2021.2.3) and commercial applications (Tobii Eye-Tracker Manager (v. 2.6.0), BrainVision Recorder (v. 1.20.0801), and BrainVision Analyzer (v. 2.3.0.8300)) to perform the experiment. We deposited the scripts and content in a GitHub Repository https://github.com/brains-on-code/AoC-FRP-Code and they are available for unrestricted open access under a CC-BY license.

Reproduced under the paper's license (CC BY), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 8 authors, 10 keywords, 6 MeSH terms, 1 funder, 32 references.

Cite

This paper

Bergum, A., Maurer, A.-M., Peitek, N., Bader, R., Mecklinger, A., Demberg, V., Siegmund, J., & Apel, S. (2026). Fixation-related potentials reveal that confusing program code elicits a late frontal positivity. Scientific reports, 16(1), 16833. https://doi.org/10.1038/s41598-026-50946-9

BibTeX

@article{bergum2026fixation,
author = {Bergum, Annabelle and Maurer, Anna-Maria and Peitek, Norman and Bader, Regine and Mecklinger, Axel and Demberg, Vera and Siegmund, Janet and Apel, Sven},
title = {{Fixation-related potentials reveal that confusing program code elicits a late frontal positivity}},
journal = {Scientific reports},
year = {2026},
month = jun,
volume = {16},
number = {1},
pages = {16833},
publisher = {Nature Publishing Group},
issn = {2045-2322},
doi = {10.1038/s41598-026-50946-9},
url = {https://doi.org/10.1038/s41598-026-50946-9},
pmid = {42225689},
pmcid = {PMC13226682}
}

RIS

TY - JOUR
AU - Bergum, Annabelle
AU - Maurer, Anna-Maria
AU - Peitek, Norman
AU - Bader, Regine
AU - Mecklinger, Axel
AU - Demberg, Vera
AU - Siegmund, Janet
AU - Apel, Sven
TI - Fixation-related potentials reveal that confusing program code elicits a late frontal positivity
T2 - Scientific reports
J2 - Sci Rep
PY - 2026
DA - 2026/06/01
VL - 16
IS - 1
SP - 16833
SN - 2045-2322
PB - Nature Publishing Group
DO - 10.1038/s41598-026-50946-9
UR - https://doi.org/10.1038/s41598-026-50946-9
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s41598-026-50946-9",
"type": "article-journal",
"title": "Fixation-related potentials reveal that confusing program code elicits a late frontal positivity",
"container-title": "Scientific reports",
"author": [
{
"family": "Bergum",
"given": "Annabelle"
},
{
"family": "Maurer",
"given": "Anna-Maria"
},
{
"family": "Peitek",
"given": "Norman"
},
{
"family": "Bader",
"given": "Regine"
},
{
"family": "Mecklinger",
"given": "Axel"
},
{
"family": "Demberg",
"given": "Vera"
},
{
"family": "Siegmund",
"given": "Janet"
},
{
"family": "Apel",
"given": "Sven"
}
],
"container-title-short": "Sci Rep",
"volume": "16",
"issue": "1",
"page": "16833",
"DOI": "10.1038/s41598-026-50946-9",
"PMID": "42225689",
"PMCID": "PMC13226682",
"ISSN": "2045-2322",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s41598-026-50946-9",
"language": "en",
"issued": {
"date-parts": [
[
2026,
6,
1
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1016/j.xcrm.2026.102766 [code]
A longitudinal single-cell and spatial multiomic atlas of pediatric high-grade glioma.
Journal: Cell reports. Medicine
In common: rpy2, lmerTest, lme4, 9 other tools
[2] doi:10.34133/csbj.0042 [code]
Using Steady-State Visual Evoked Potentials to Characterize Wide-Ranging Retinopathy Linked to &lt;i&gt;CRB1&lt;/i&gt;: Implications for Clinical Trials.
Journal: Computational and structural biotechnology journal
In common: PsychoPy, MNE-Python, lmerTest, 8 other tools, EEG
[3] doi:10.1038/s41598-026-41129-7 [code]
Non-linear relationships between auditory mismatch responses and the inharmonicity of complex sounds.
Journal: Scientific reports
In common: PsychoPy, nlme, MNE-Python, 7 other tools, EEG
[4] doi:10.1038/s41597-025-05174-7 [code]
A large-scale MEG and EEG dataset for object recognition in naturalistic scenes
Journal: n/a
In common: MNE-Python, OpenCV, h5py, 7 other tools, methods / tools, EEG
[5] doi:10.21203/rs.3.rs-9676637/v1 [code]
A Comprehensive Benchmarking of Spatial Deconvolution and Domain Detection Methods across Diverse Tissues and Spatial Transcriptomic Technologies
Journal: Research Square (preprint)
In common: rpy2, OpenCV, h5py, 7 other tools, methods / tools
[6] doi:10.1038/s42003-026-10957-8 [code]
Brain defence by the extracellular matrix protein Cochlin.
Journal: Communications biology
In common: nlme, lme4, OpenCV, 7 other tools
[7] doi:10.7554/elife.108673 [code]
Adaptive behavior is guided by integrated representations of controlled and non-controlled information.
Journal: eLife
In common: rpy2, MNE-Python, lme4, 6 other tools, EEG
[8] doi:10.1117/1.nph.13.2.025001 [code]
Surface-based image reconstruction optimization for high-density functional near-infrared spectroscopy.
Journal: Neurophotonics
In common: MNE-Python, OpenCV, h5py, 7 other tools
[9] doi:10.1093/cercor/bhag113 [code]
Long-term reliability and stability of parameterized resting state EEG: evidence from a five-year follow-up.
Journal: Cerebral cortex (New York, N.Y. : 1991)
In common: MNE-Python, lmerTest, lme4, 6 other tools, methods / tools, EEG
[10] doi:10.1371/journal.pone.0343722 [code]
Comprehensive methodology for sample enrichment in EEG biomarker studies for Alzheimer's risk classification.
Journal: PloS one
In common: rpy2, MNE-Python, Pillow, 6 other tools, EEG

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.