feat: add require_matched_subjects flag and improve difference wave reporting

- Add explicit control for difference wave subject matching
- Update write_processing_report to display per-code thresholds
- Track subjects contributing to each condition
- Pass min_trials_per_code to report generation

Author: Buch

code/postprocessing/postprocessing-eeg/batch_eeg_postprocessing.m

@@ -74,8 +74,8 @@
 
 
 % inclusion thresholds
-min_epochs_threshold = 10;
 min_accuracy_threshold = 0.6;
+min_epochs_threshold = 10;
 
 % RT trimming parameters
 rt_lower_bound = 150; % ms, set to 0 if no lower bound
@@ -119,6 +119,9 @@
 
 fprintf('defined %d difference wave computations\n', size(diff_waves_table, 1));
 
+% difference wave computation settings
+require_matched_subjects = true;  % only compute difference waves from subjects present in both conditions
+
 ern_time_window = [0, 100];
 pe_time_window = [200, 500];
 
@@ -139,7 +142,6 @@
 processing_report.input_data.code_names = code_names;
 
 % processing parameters
-processing_report.parameters.min_epochs_threshold = min_epochs_threshold;
 processing_report.parameters.min_accuracy_threshold = min_accuracy_threshold;
 processing_report.parameters.rt_lower_bound = rt_lower_bound;
 processing_report.parameters.rt_outlier_threshold = rt_outlier_threshold;
@@ -190,7 +192,7 @@
     % generate subject summary table
     fprintf('\ngenerating subject summary table...\n');
     grand_avg_data = load(grand_avg_file);
-    generate_subject_summary_table(grand_avg_data.processing_stats, included_subjects, codes, output_dir, grand_avg_data.condition_inclusion);
+    generate_subject_summary_table(grand_avg_data.processing_stats, included_subjects, codes, output_dir, grand_avg_data.condition_inclusion, grand_avg_data.min_trials_per_code);
 
 else
     % load existing grand averages if not running step 1
@@ -217,7 +219,7 @@
 
     % function call
     grand_avg_data = load(grand_avg_file);
-    [difference_waves] = compute_difference_waves(grand_averages, included_subjects, diff_waves_table, output_dir, grand_avg_data.condition_inclusion, codes);
+    [difference_waves] = compute_difference_waves(grand_averages, included_subjects, diff_waves_table, output_dir, grand_avg_data.condition_inclusion, codes, require_matched_subjects);
 
     step2_end_time = datetime('now');
 
@@ -278,8 +280,18 @@
 
 % generate comprehensive text report
 report_file = fullfile(output_dir, sprintf('postprocessing_summary_%s.txt', datestr(now, 'yyyy-mm-dd_HH-MM-SS')));
-write_processing_report(processing_report, diff_waves_table, output_dir, report_file);
-
+% load min_trials_per_code if available
+if exist(grand_avg_file, 'file')
+    grand_avg_data = load(grand_avg_file);
+    if isfield(grand_avg_data, 'min_trials_per_code')
+        min_trials_per_code = grand_avg_data.min_trials_per_code;
+    else
+        min_trials_per_code = [];
+    end
+else
+    min_trials_per_code = [];
+end
+write_processing_report(processing_report, diff_waves_table, output_dir, report_file, min_trials_per_code);
 % print console summary
 print_console_summary(processing_report, output_dir, report_file, console_log_file);
 

code/postprocessing/postprocessing-eeg/compute_difference_waves.m

@@ -1,4 +1,4 @@
-function difference_waves = compute_difference_waves(grand_averages, included_subjects, diff_waves_table, output_dir, condition_inclusion, codes)
+function difference_waves = compute_difference_waves(grand_averages, included_subjects, diff_waves_table, output_dir, condition_inclusion, codes, require_matched_subjects)
 % compute_difference_waves - create configurable difference waves for ERP analysis
 %
 % this function computes difference waves based on user-defined table
@@ -8,6 +8,7 @@
 %   included_subjects - cell array of included subject IDs
 %   diff_waves_table - matrix [minuend_code, subtrahend_code, wave_name]
 %   output_dir - path to save difference waves
+%   require_matched_subjects - logical, if true only compute from subjects in both conditions (default: true)
 %
 % outputs:
 %   difference_waves - struct containing computed difference waves
@@ -59,6 +60,7 @@
 fprintf(log_fid, 'script: compute_difference_waves.m\n');
 fprintf(log_fid, 'included subjects: %d\n', length(included_subjects));
 fprintf(log_fid, 'subjects: %s\n', strjoin(included_subjects, ', '));
+fprintf(log_fid, 'require_matched_subjects: %s\n', mat2str(require_matched_subjects));
 fprintf(log_fid, '\n');
 
 % initialize difference waves structure
@@ -107,30 +109,36 @@
             subject = included_subjects{subj_idx};
             subject_inclusion = condition_inclusion(subject);
 
-            % check if subject has both conditions
-            if subject_inclusion(minuend_idx) && subject_inclusion(subtrahend_idx)
-                subjects_for_this_diffwave{end+1} = subject;
-
-                % find this subject's position in each grand average
-                % count how many subjects before this one have the minuend code
-                minuend_pos = 0;
-                for i = 1:subj_idx
-                    temp = condition_inclusion(included_subjects{i});
-                    if temp(minuend_idx)
-                        minuend_pos = minuend_pos + 1;
+            % check if subject has both conditions (only if require_matched_subjects is true)
+            if require_matched_subjects
+                % strict mode: subject must have BOTH conditions
+                if subject_inclusion(minuend_idx) && subject_inclusion(subtrahend_idx)
+                    subjects_for_this_diffwave{end+1} = subject;
+
+                    % find this subject's position in each grand average
+                    % count how many subjects before this one have the minuend code
+                    minuend_pos = 0;
+                    for i = 1:subj_idx
+                        temp = condition_inclusion(included_subjects{i});
+                        if temp(minuend_idx)
+                            minuend_pos = minuend_pos + 1;
+                        end
                     end
-                end
 
-                subtrahend_pos = 0;
-                for i = 1:subj_idx
-                    temp = condition_inclusion(included_subjects{i});
-                    if temp(subtrahend_idx)
-                        subtrahend_pos = subtrahend_pos + 1;
+                    subtrahend_pos = 0;
+                    for i = 1:subj_idx
+                        temp = condition_inclusion(included_subjects{i});
+                        if temp(subtrahend_idx)
+                            subtrahend_pos = subtrahend_pos + 1;
+                        end
                     end
-                end
 
-                minuend_indices = [minuend_indices, minuend_pos];
-                subtrahend_indices = [subtrahend_indices, subtrahend_pos];
+                    minuend_indices = [minuend_indices, minuend_pos];
+                    subtrahend_indices = [subtrahend_indices, subtrahend_pos];
+                end
+            else
+                % loose mode would go here if ever needed
+                error('Loose mode (require_matched_subjects=false) not implemented for within-subjects design');
             end
         end
 

code/postprocessing/postprocessing-eeg/generate_subject_summary_table.m

@@ -1,4 +1,4 @@
-function generate_subject_summary_table(processing_stats, included_subjects, codes, output_dir, condition_inclusion)
+function generate_subject_summary_table(processing_stats, included_subjects, codes, output_dir, condition_inclusion, min_trials_per_code)
 % generate tab-separated summary table for all subjects
 % creates table matching format from google doc for easy copy-paste
 %

code/postprocessing/postprocessing-eeg/write_processing_report.m

@@ -1,4 +1,4 @@
-function write_processing_report(processing_report, diff_waves_table, output_dir, report_file)
+function write_processing_report(processing_report, diff_waves_table, output_dir, report_file, min_trials_per_code)
 % write comprehensive processing summary report to text file
 % 
 % inputs:
@@ -43,7 +43,15 @@ function write_processing_report(processing_report, diff_waves_table, output_dir
 
 % processing parameters
 fprintf(report_fid, '=== PROCESSING PARAMETERS ===\n');
-fprintf(report_fid, 'minimum epochs threshold: %d\n', processing_report.parameters.min_epochs_threshold);
+fprintf(report_fid, 'minimum epochs thresholds (per code):\n');
+if ~isempty(min_trials_per_code)
+    codes_list = fieldnames(min_trials_per_code);
+    for i = 1:length(codes_list)
+        fprintf(report_fid, '  - %s: %d trials\n', codes_list{i}, min_trials_per_code.(codes_list{i}));
+    end
+else
+    fprintf(report_fid, '  - not available\n');
+end
 fprintf(report_fid, 'minimum accuracy threshold: %.1f%%\n', processing_report.parameters.min_accuracy_threshold * 100);
 fprintf(report_fid, 'RT lower bound: %d ms\n', processing_report.parameters.rt_lower_bound);
 fprintf(report_fid, 'RT outlier threshold: %.1f SD\n', processing_report.parameters.rt_outlier_threshold);