Merge pull request #1 from barbaLab/copilot/add-modulation-index-script

Add d-prime modulation index metric to +metrics package

Author: Nabarb

+metrics/+compute/dprime.m

@@ -0,0 +1,49 @@
+function dp = dprime(data, baseline_idx, signal_idx)
+%DPRIME computes the d-prime modulation index between two time windows
+%
+% Inputs:
+%   data          - nDims x nTime x nTrials array of embedded signals
+%   baseline_idx  - indices of baseline time points
+%   signal_idx    - indices of signal time points
+%
+% Output:
+%   dp            - 1 x nDims vector of d' values
+
+if isempty(data)
+    dp = nan;
+    return;
+end
+
+[nDims, ~, nTrials] = size(data);
+
+if nTrials < 2
+    dp = nan(1, nDims);
+    return;
+end
+
+% Mean activity per trial within each window
+baseline = squeeze(mean(data(:, baseline_idx, :), 2));  % nDims x nTrials
+signal   = squeeze(mean(data(:, signal_idx,   :), 2));  % nDims x nTrials
+
+% When nDims == 1, squeeze reduces 1x1xnTrials to a nTrials×1 column;
+% reshape to 1×nTrials so that trailing dimension is always trials.
+if nDims == 1
+    baseline = baseline(:)';
+    signal   = signal(:)';
+end
+
+% Compute d' per dimension (statistics across trials)
+mu_b  = mean(baseline, 2);
+mu_s  = mean(signal,   2);
+std_b = std(baseline,  0, 2);
+std_s = std(signal,    0, 2);
+
+% Pooled standard deviation: sqrt( (sigma_b^2 + sigma_s^2) / 2 )
+pooled_std = sqrt(0.5 .* (std_b.^2 + std_s.^2));
+
+dp = (mu_s - mu_b) ./ pooled_std;
+
+% Return as 1 x nDims row vector
+dp = dp';
+
+end

+metrics/+pars/DPrime.m

@@ -0,0 +1,7 @@
+function pars = DPrime()
+
+pars.baseline_idx = [];
+pars.signal_idx   = [];
+
+
+end

+metrics/DPrime.m

@@ -0,0 +1,18 @@
+function dp = DPrime(E, pars)
+%DPRIME computes d-prime modulation index between baseline and signal windows
+%
+% input
+%                       E             1xnTrials cell array, each with nDims×nTime matrix
+%                       pars          parameter struct with fields:
+%                                       .baseline_idx  indices of baseline time points
+%                                       .signal_idx    indices of signal time points
+% output
+%                       dp            1×nDims vector of d' values
+
+if nargin < 2 || isempty(pars)
+    pars = metrics.pars.DPrime();
+end
+
+dp = metrics.compute.dprime(cat(3, E{:}), pars.baseline_idx, pars.signal_idx);
+
+end

@NeuralEmbedding/NeuralEmbedding.m

@@ -69,6 +69,8 @@
 
         % Metrics computation
         appendM                 = false
+        baseline_idx            = []
+        signal_idx              = []
 
         % Metadata
         Meta                    = struct("AnName","","ExpGroup","");

@NeuralEmbedding/computeMetrics.m

@@ -19,6 +19,10 @@
         case {'Alignment','alignment','align','Align'}
             type = "Alignment";
             parNames = [""];
+
+        case {'DPrime','dprime','dPrime','d_prime','d-prime'}
+            type = "DPrime";
+            parNames = ["baseline_idx","signal_idx"];
         otherwise
             error("Specific metric not found.")
     end