final changes for new release

This commit is contained in:
2026-03-28 14:30:20 -07:00
parent 70c4c1e779
commit 74ce2eeb2e
3 changed files with 158 additions and 67 deletions
+96 -45
View File
@@ -149,9 +149,23 @@ PSP: bool
PSP_TIME_WINDOW: int
PSP_THRESHOLD: float
CV: bool
CV_THRESHOLD: int
MAD: bool
MAD_THRESHOLD: int
PSD_NOISE: bool
TARGET_FREQ_DIV: int
DB_LIMIT: int
CHANNEL_VAR: bool
CHANNEL_THRESH: float
BAD_CHANNELS_HANDLING: str
MAX_DIST: float
MIN_NEIGHBORS: int
MAX_BAD_CHANNELS: int
TDDR: bool
@@ -2974,7 +2988,20 @@ def run_second_level_analysis(df_contrasts, raw, p, bounds):
'mean_beta': mean_beta,
'n_subjects': len(Y)
})
if not group_results:
# Create a "Warning" figure instead of a map
fig, ax = plt.subplots(figsize=(8, 4))
ax.text(0.5, 0.5,
f"Second-Level Analysis Aborted\n\n"
f"Reason: All {len(channels)} channels skipped.\n"
f"Requirement: At least 2 subjects (IDs) per channel.\n"
f"Current Subject Count: {df_contrasts['ID'].nunique()}",
ha='center', va='center', fontsize=12, color='darkred',
bbox=dict(facecolor='white', alpha=0.5, edgecolor='red'))
ax.set_axis_off()
plt.show()
return pd.DataFrame()
df_group = pd.DataFrame(group_results)
logger.info("Second-level results:\n%s", df_group)
@@ -3652,14 +3679,14 @@ def detect_sensor_displacement(raw, threshold_ratio=0.05):
def detect_high_freq_noise(raw, db_limit=-60):
def detect_high_freq_noise(raw, db_limit=-60, freq_div=4):
"""
Identifies channels with excessive power at high frequencies
(sfreq/4), usually indicating electronic interference.
"""
ch_names = raw.ch_names
sfreq = raw.info['sfreq']
target_freq = sfreq / 4
target_freq = sfreq / freq_div
# Compute PSD
spectrum = raw.compute_psd(fmin=0.1, fmax=sfreq/2)
@@ -3952,18 +3979,33 @@ def process_participant(file_path, progress_callback=None):
if progress_callback: progress_callback(8)
logger.info("Step 8 Completed.")
# TODO: Add callbacks and user defined parameters
bad_cv, fig_cv = find_bad_channels_cv(raw, cv_threshold=20.0)
fig_individual['cv'] = fig_cv
bad_cv = []
if CV:
bad_cv, fig_cv = find_bad_channels_cv(raw, cv_threshold=CV_THRESHOLD)
fig_individual['cv'] = fig_cv
if progress_callback: progress_callback(9)
logger.info("Step 9 Completed.")
bad_range, fig_range = find_bad_channels_range(raw, threshold=3.0)
fig_individual['range'] = fig_range
bad_range = []
if MAD:
bad_range, fig_range = find_bad_channels_range(raw, threshold=MAD_THRESHOLD)
fig_individual['range'] = fig_range
if progress_callback: progress_callback(10)
logger.info("Step 10 Completed.")
bad_noise, fig_noise = detect_high_freq_noise(raw, db_limit=-60)
fig_individual['psd_noise'] = fig_noise
bad_noise = []
if PSD_NOISE:
bad_noise, fig_noise = detect_high_freq_noise(raw, db_limit=DB_LIMIT, freq_div=TARGET_FREQ_DIV)
fig_individual['psd_noise'] = fig_noise
if progress_callback: progress_callback(11)
logger.info("Step 11 Completed.")
bad_disp, fig_disp = detect_sensor_displacement(raw, threshold_ratio=0.05)
fig_individual['displacement'] = fig_disp
bad_disp = []
if CHANNEL_VAR:
bad_disp, fig_disp = detect_sensor_displacement(raw, threshold_ratio=CHANNEL_THRESH)
fig_individual['displacement'] = fig_disp
if progress_callback: progress_callback(12)
logger.info("Step 12 Completed.")
# Step 9: Bad Channels Handling
if BAD_CHANNELS_HANDLING != "None":
@@ -3977,77 +4019,86 @@ def process_participant(file_path, progress_callback=None):
fig_individual["fig4"] = fig_raw_after
fig_individual["Compare"] = fig_compare
elif BAD_CHANNELS_HANDLING == "Remove":
#NOTE: testing this
num_bad = len(bad_channels)
# Check against the threshold
if num_bad > MAX_BAD_CHANNELS:
raise Exception(
f"Data Quality Error: {num_bad} channels flagged for removal, "
f"which exceeds the limit of {MAX_BAD_CHANNELS}. To avoid this,"
f"either lower your filtering parameters or increase MAX_BAD_CHANNELS."
)
raw.pick_types(fnirs=True, exclude='bads')
logger.info(f"Physically removed {len(bad_channels)} channels from the dataset.")
if progress_callback: progress_callback(9)
logger.info("Step 9 Completed.")
if progress_callback: progress_callback(13)
logger.info("Step 13 Completed.")
# Step 10: Optical Density
raw_od = optical_density(raw)
fig_raw_od = raw_od.plot(duration=raw.times[-1], n_channels=raw.info['nchan'], title="Optical Density", show=False)
fig_individual["Optical Density"] = fig_raw_od
if progress_callback: progress_callback(10)
logger.info("Step 10 Completed.")
if progress_callback: progress_callback(14)
logger.info("Step 14 Completed.")
# Step 11: Temporal Derivative Distribution Repair Filtering
if TDDR:
raw_od = temporal_derivative_distribution_repair(raw_od)
fig_raw_od_tddr = raw_od.plot(duration=raw.times[-1], n_channels=raw.info['nchan'], title="After TDDR (Motion Correction)", show=False)
fig_individual["TDDR"] = fig_raw_od_tddr
if progress_callback: progress_callback(11)
logger.info("Step 11 Completed.")
if progress_callback: progress_callback(15)
logger.info("Step 15 Completed.")
# Step 12: Wavelet Filtering
if WAVELET:
raw_od, fig = calculate_and_apply_wavelet(raw_od)
fig_individual["Wavelet"] = fig
if progress_callback: progress_callback(12)
logger.info("Step 12 Completed.")
if progress_callback: progress_callback(16)
logger.info("Step 16 Completed.")
# Step 13: Haemoglobin Concentration
raw_haemo = beer_lambert_law(raw_od, ppf=calculate_dpf(file_path))
fig_raw_haemo_bll = raw_haemo.plot(duration=raw_haemo.times[-1], n_channels=raw_haemo.info['nchan'], title="HbO and HbR Signals", show=False)
fig_individual["BLL"] = fig_raw_haemo_bll
if progress_callback: progress_callback(13)
logger.info("Step 13 Completed.")
if progress_callback: progress_callback(17)
logger.info("Step 17 Completed.")
# Step 14: Enhance Negative Correlation
if ENHANCE_NEGATIVE_CORRELATION:
raw_haemo = enhance_negative_correlation(raw_haemo)
fig_raw_haemo_enc = raw_haemo.plot(duration=raw_haemo.times[-1], n_channels=raw_haemo.info['nchan'], title="Enhance Negative Correlation", show=False)
fig_individual["ENC"] = fig_raw_haemo_enc
if progress_callback: progress_callback(14)
logger.info("Step 14 Completed.")
if progress_callback: progress_callback(18)
logger.info("Step 18 Completed.")
# Step 15: Filter
if FILTER:
raw_haemo, fig_filter, fig_raw_haemo_filter = filter_the_data(raw_haemo)
fig_individual["filter1"] = fig_filter
fig_individual["filter2"] = fig_raw_haemo_filter
if progress_callback: progress_callback(15)
logger.info("Step 15 Completed.")
if progress_callback: progress_callback(19)
logger.info("Step 19 Completed.")
# Step 16: Extracting Events
events, event_dict = events_from_annotations(raw_haemo)
fig_events = plot_events(events, event_id=event_dict, sfreq=raw_haemo.info["sfreq"], show=False)
fig_individual["events"] = fig_events
if progress_callback: progress_callback(16)
logger.info("Step 16 Completed.")
if progress_callback: progress_callback(20)
logger.info("Step 20 Completed.")
# Step 17: Epoch Calculations
epochs, fig_epochs = epochs_calculations(raw_haemo, events, event_dict)
for name, fig in fig_epochs:
fig_individual[f"epochs_{name}"] = fig
if progress_callback: progress_callback(17)
logger.info("Step 17 Completed.")
if progress_callback: progress_callback(21)
logger.info("Step 21 Completed.")
# Step 18: Design Matrix
design_matrix, fig_design_matrix = make_design_matrix(raw_haemo, short_chans)
fig_individual["Design Matrix"] = fig_design_matrix
if progress_callback: progress_callback(18)
logger.info("Step 18 Completed.")
if progress_callback: progress_callback(22)
logger.info("Step 22 Completed.")
# Step 19: Run GLM
@@ -4063,24 +4114,24 @@ def process_participant(file_path, progress_callback=None):
# A large p-value means the data do not provide strong evidence that the effect is different from zero.
if progress_callback: progress_callback(19)
logger.info("19")
if progress_callback: progress_callback(23)
logger.info("23")
# Step 20: Generate GLM Results
if "derivative" not in HRF_MODEL.lower():
fig_glm_result = plot_glm_results(file_path, raw_haemo, glm_est, design_matrix)
for name, fig in fig_glm_result:
fig_individual[f"GLM {name}"] = fig
if progress_callback: progress_callback(20)
logger.info("20")
if progress_callback: progress_callback(24)
logger.info("24")
# Step 21: Generate Channel Significance
if HRF_MODEL == "fir":
fig_significance = individual_significance(raw_haemo, glm_est)
for name, fig in fig_significance:
fig_individual[f"Significance {name}"] = fig
if progress_callback: progress_callback(21)
logger.info("21")
if progress_callback: progress_callback(25)
logger.info("25")
# Step 22: Generate Channel, Region of Interest, and Contrast Results
cha = glm_est.to_dataframe()
@@ -4136,8 +4187,8 @@ def process_participant(file_path, progress_callback=None):
contrast_dict[condition] = contrast_vector
if progress_callback: progress_callback(22)
logger.info("22")
if progress_callback: progress_callback(26)
logger.info("26")
# Step 23: Compute Contrast Results
contrast_results = {}
@@ -4152,16 +4203,16 @@ def process_participant(file_path, progress_callback=None):
cha["ID"] = file_path
if progress_callback: progress_callback(23)
logger.info("23")
if progress_callback: progress_callback(27)
logger.info("27")
# Step 24: Finishing Up
fig_bytes = convert_fig_dict_to_png_bytes(fig_individual)
sanitize_paths_for_pickle(raw_haemo, epochs)
if progress_callback: progress_callback(25)
logger.info("25")
if progress_callback: progress_callback(28)
logger.info("28")
# TODO: Tidy up
# Extract the parameters this file was ran with. No need to return age, gender, group?