From 1af1a723fc7fdd4f51073a7b974ea6011b971b99 Mon Sep 17 00:00:00 2001 From: Micha Burkhardt Date: Fri, 14 Feb 2025 10:55:29 +0100 Subject: [PATCH 1/8] added std parameter to SW methods and limit the width to the window size --- pydfc/dfc_methods/discrete_hmm.py | 1 + pydfc/dfc_methods/sliding_window.py | 35 ++++++++-------------- pydfc/dfc_methods/sliding_window_clustr.py | 1 + 3 files changed, 15 insertions(+), 22 deletions(-) diff --git a/pydfc/dfc_methods/discrete_hmm.py b/pydfc/dfc_methods/discrete_hmm.py index 0bb6779..e57de21 100644 --- a/pydfc/dfc_methods/discrete_hmm.py +++ b/pydfc/dfc_methods/discrete_hmm.py @@ -65,6 +65,7 @@ def __init__(self, **params): "backend", "n_subj_clstrs", "W", + "window_std", "n_overlap", "n_states", "normalization", diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 17ec5e6..0fcedbf 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -45,6 +45,7 @@ def __init__(self, **params): "is_state_based", "sw_method", "tapered_window", + "window_std", "W", "n_overlap", "normalization", @@ -127,7 +128,7 @@ def FC(self, time_series): return C - def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False): + def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False, window_std=None): # W is in time samples L = time_series.shape[1] @@ -141,36 +142,25 @@ def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False): model.fit(time_series.T) self.graphical_lasso_alpha_ = model.alpha_ - window_taper = signal.windows.gaussian(W, std=3 * W / 22) - # C = DFC(measure=self) FCSs = list() TR_array = list() for l in range(0, L - W + 1, step): - ######### creating a rectangel window ############ + # Create rectangular window window = np.zeros((L)) window[l : l + W] = 1 - ########### tapering the window ############## + # Taper the window if tapered_window: - window = signal.convolve(window, window_taper, mode="same") / sum( - window_taper - ) - - window = np.repeat( - np.expand_dims(window, axis=0), time_series.shape[0], axis=0 - ) - - # int(l-W/2):int(l+3*W/2) is the nonzero interval after tapering - FCSs.append( - self.FC( - np.multiply(time_series, window)[ - :, max(int(l - W / 2), 0) : min(int(l + 3 * W / 2), L) - ] - ) - ) + std = window_std if window_std is not None else 3 * W / 22 + window_taper = signal.windows.gaussian(W, std=std) + #window = signal.convolve(window, window_taper, mode="same") / sum(window_taper) + window[l:l+W] = window_taper + + window = np.repeat(np.expand_dims(window, axis=0), time_series.shape[0], axis=0) + FCSs.append(self.FC(np.multiply(time_series, window)[:,l:l+W])) TR_array.append(int((l + (l + W)) / 2)) - + return np.array(FCSs), np.array(TR_array) def estimate_FCS(self, time_series): @@ -200,6 +190,7 @@ def estimate_dFC(self, time_series): W=int(self.params["W"] * time_series.Fs), n_overlap=self.params["n_overlap"], tapered_window=self.params["tapered_window"], + window_std=self.params["window_std"] ) # record time diff --git a/pydfc/dfc_methods/sliding_window_clustr.py b/pydfc/dfc_methods/sliding_window_clustr.py index c19e883..137552c 100644 --- a/pydfc/dfc_methods/sliding_window_clustr.py +++ b/pydfc/dfc_methods/sliding_window_clustr.py @@ -67,6 +67,7 @@ def __init__(self, clstr_distance="euclidean", **params): "coi_correction", "n_subj_clstrs", "W", + "window_std", "n_overlap", "n_states", "normalization", From fd8db391bba1cd977f220e57aedf52fc80b08dc9 Mon Sep 17 00:00:00 2001 From: Micha Burkhardt Date: Fri, 14 Feb 2025 10:56:10 +0100 Subject: [PATCH 2/8] vectorized pearson correlation --- pydfc/dfc_methods/sliding_window.py | 37 ++++++++++------------------- 1 file changed, 13 insertions(+), 24 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 0fcedbf..8ab5e2d 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -96,35 +96,24 @@ def calc_MI(self, X, Y): return MI def FC(self, time_series): - - if self.params["sw_method"] == "GraphLasso": - # Standardize the data (zero mean, unit variance for each feature) - mean = np.mean(time_series, axis=1, keepdims=True) - std = np.std(time_series, axis=1, keepdims=True) - time_series_standardized = np.where(std != 0, (time_series - mean) / std, 0) - model = GraphicalLasso(alpha=self.graphical_lasso_alpha_) - model.fit(time_series_standardized.T) - # the covariance matrix will equal the correlation matrix + # Graphical Lasso + if self.params['sw_method']=='GraphLasso': + model = GraphicalLassoCV() + model.fit(time_series.T) C = model.covariance_ - else: + # Mutual information + elif self.params['sw_method']=='MI': C = np.zeros((time_series.shape[0], time_series.shape[0])) + for i in range(time_series.shape[0]): - for j in range(i, time_series.shape[0]): - + for j in range(i, time_series.shape[0]): X = time_series[i, :] Y = time_series[j, :] - - if self.params["sw_method"] == "MI": - ########### Mutual Information ############## - C[j, i] = self.calc_MI(X, Y) - else: - ########### Pearson Correlation ############## - if np.var(X) == 0 or np.var(Y) == 0: - C[j, i] = 0 - else: - C[j, i] = np.corrcoef(X, Y)[0, 1] - - C[i, j] = C[j, i] + C[j, i] = self.calc_MI(X, Y) + # Pearson correlation + else: + C = np.corrcoef(time_series) + C[np.isnan(C)] = 0 return C From cd270e35981f393bd5ff7aa24a3420cbfeea1af8 Mon Sep 17 00:00:00 2001 From: Micha Burkhardt Date: Tue, 18 Feb 2025 13:53:45 +0100 Subject: [PATCH 3/8] Set diagonal to 1 Co-authored-by: MOHAMMAD TORABI <72619172+mtorabi59@users.noreply.github.com> --- pydfc/dfc_methods/sliding_window.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 8ab5e2d..05266a5 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -114,7 +114,8 @@ def FC(self, time_series): else: C = np.corrcoef(time_series) C[np.isnan(C)] = 0 - + # make the diagonal elements 1 (for nan values on the diagonal) + C[np.diag_indices_from(C)] = 1 return C def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False, window_std=None): From 0702f257b537801a1a6f36e84d04dedb20425ec9 Mon Sep 17 00:00:00 2001 From: Micha Burkhardt Date: Tue, 18 Feb 2025 13:55:16 +0100 Subject: [PATCH 4/8] readded accidentally removed lines for GraphLasso --- pydfc/dfc_methods/sliding_window.py | 11 +++++++++-- 1 file changed, 9 insertions(+), 2 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 05266a5..2f6620c 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -98,9 +98,15 @@ def calc_MI(self, X, Y): def FC(self, time_series): # Graphical Lasso if self.params['sw_method']=='GraphLasso': - model = GraphicalLassoCV() - model.fit(time_series.T) + # Standardize the data (zero mean, unit variance for each feature) + mean = np.mean(time_series, axis=1, keepdims=True) + std = np.std(time_series, axis=1, keepdims=True) + time_series_standardized = np.where(std != 0, (time_series - mean) / std, 0) + model = GraphicalLasso(alpha=self.graphical_lasso_alpha_) + model.fit(time_series_standardized.T) + # the covariance matrix will equal the correlation matrix C = model.covariance_ + # Mutual information elif self.params['sw_method']=='MI': C = np.zeros((time_series.shape[0], time_series.shape[0])) @@ -110,6 +116,7 @@ def FC(self, time_series): X = time_series[i, :] Y = time_series[j, :] C[j, i] = self.calc_MI(X, Y) + # Pearson correlation else: C = np.corrcoef(time_series) From 40d2031cddedbc74aa1c5bc2f20f0d31be36f73b Mon Sep 17 00:00:00 2001 From: Micha Burkhardt Date: Tue, 18 Feb 2025 14:05:42 +0100 Subject: [PATCH 5/8] readded convolution to follow Allen et al. --- pydfc/dfc_methods/sliding_window.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 2f6620c..84a2918 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -151,8 +151,7 @@ def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False, window_ if tapered_window: std = window_std if window_std is not None else 3 * W / 22 window_taper = signal.windows.gaussian(W, std=std) - #window = signal.convolve(window, window_taper, mode="same") / sum(window_taper) - window[l:l+W] = window_taper + window = signal.convolve(window, window_taper, mode="same") / sum(window_taper) window = np.repeat(np.expand_dims(window, axis=0), time_series.shape[0], axis=0) FCSs.append(self.FC(np.multiply(time_series, window)[:,l:l+W])) From 9f10bf7145560e1387d83f96cd7aa1123d815527 Mon Sep 17 00:00:00 2001 From: MOHAMMAD TORABI <72619172+mtorabi59@users.noreply.github.com> Date: Tue, 4 Mar 2025 14:54:46 -0500 Subject: [PATCH 6/8] minor changes --- pydfc/dfc_methods/sliding_window.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 84a2918..7ad7e68 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -97,7 +97,7 @@ def calc_MI(self, X, Y): def FC(self, time_series): # Graphical Lasso - if self.params['sw_method']=='GraphLasso': + if self.params["sw_method"]=="GraphLasso": # Standardize the data (zero mean, unit variance for each feature) mean = np.mean(time_series, axis=1, keepdims=True) std = np.std(time_series, axis=1, keepdims=True) @@ -108,7 +108,7 @@ def FC(self, time_series): C = model.covariance_ # Mutual information - elif self.params['sw_method']=='MI': + elif self.params["sw_method"]=="MI": C = np.zeros((time_series.shape[0], time_series.shape[0])) for i in range(time_series.shape[0]): From c3b7579a663e5046e8fdbed78113fc3b657f7047 Mon Sep 17 00:00:00 2001 From: MOHAMMAD TORABI <72619172+mtorabi59@users.noreply.github.com> Date: Tue, 4 Mar 2025 15:03:06 -0500 Subject: [PATCH 7/8] minor changes --- pydfc/dfc_methods/sliding_window.py | 22 ++++++++++++++-------- 1 file changed, 14 insertions(+), 8 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 7ad7e68..0d229ef 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -97,7 +97,7 @@ def calc_MI(self, X, Y): def FC(self, time_series): # Graphical Lasso - if self.params["sw_method"]=="GraphLasso": + if self.params["sw_method"] == "GraphLasso": # Standardize the data (zero mean, unit variance for each feature) mean = np.mean(time_series, axis=1, keepdims=True) std = np.std(time_series, axis=1, keepdims=True) @@ -108,11 +108,11 @@ def FC(self, time_series): C = model.covariance_ # Mutual information - elif self.params["sw_method"]=="MI": + elif self.params["sw_method"] == "MI": C = np.zeros((time_series.shape[0], time_series.shape[0])) for i in range(time_series.shape[0]): - for j in range(i, time_series.shape[0]): + for j in range(i, time_series.shape[0]): X = time_series[i, :] Y = time_series[j, :] C[j, i] = self.calc_MI(X, Y) @@ -125,7 +125,9 @@ def FC(self, time_series): C[np.diag_indices_from(C)] = 1 return C - def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False, window_std=None): + def dFC( + self, time_series, W=None, n_overlap=None, tapered_window=False, window_std=None + ): # W is in time samples L = time_series.shape[1] @@ -151,10 +153,14 @@ def dFC(self, time_series, W=None, n_overlap=None, tapered_window=False, window_ if tapered_window: std = window_std if window_std is not None else 3 * W / 22 window_taper = signal.windows.gaussian(W, std=std) - window = signal.convolve(window, window_taper, mode="same") / sum(window_taper) + window = signal.convolve(window, window_taper, mode="same") / sum( + window_taper + ) - window = np.repeat(np.expand_dims(window, axis=0), time_series.shape[0], axis=0) - FCSs.append(self.FC(np.multiply(time_series, window)[:,l:l+W])) + window = np.repeat( + np.expand_dims(window, axis=0), time_series.shape[0], axis=0 + ) + FCSs.append(self.FC(np.multiply(time_series, window)[:, l : l + W])) TR_array.append(int((l + (l + W)) / 2)) return np.array(FCSs), np.array(TR_array) @@ -186,7 +192,7 @@ def estimate_dFC(self, time_series): W=int(self.params["W"] * time_series.Fs), n_overlap=self.params["n_overlap"], tapered_window=self.params["tapered_window"], - window_std=self.params["window_std"] + window_std=self.params["window_std"], ) # record time From 6ba60e60e96305ddebc8f2707d9f0a0f83fbd047 Mon Sep 17 00:00:00 2001 From: MOHAMMAD TORABI <72619172+mtorabi59@users.noreply.github.com> Date: Tue, 4 Mar 2025 15:08:03 -0500 Subject: [PATCH 8/8] minor changes --- pydfc/dfc_methods/sliding_window.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/pydfc/dfc_methods/sliding_window.py b/pydfc/dfc_methods/sliding_window.py index 0d229ef..ac202fd 100644 --- a/pydfc/dfc_methods/sliding_window.py +++ b/pydfc/dfc_methods/sliding_window.py @@ -106,17 +106,17 @@ def FC(self, time_series): model.fit(time_series_standardized.T) # the covariance matrix will equal the correlation matrix C = model.covariance_ - + # Mutual information elif self.params["sw_method"] == "MI": C = np.zeros((time_series.shape[0], time_series.shape[0])) - + for i in range(time_series.shape[0]): for j in range(i, time_series.shape[0]): X = time_series[i, :] Y = time_series[j, :] C[j, i] = self.calc_MI(X, Y) - + # Pearson correlation else: C = np.corrcoef(time_series) @@ -156,13 +156,13 @@ def dFC( window = signal.convolve(window, window_taper, mode="same") / sum( window_taper ) - + window = np.repeat( np.expand_dims(window, axis=0), time_series.shape[0], axis=0 ) FCSs.append(self.FC(np.multiply(time_series, window)[:, l : l + W])) TR_array.append(int((l + (l + W)) / 2)) - + return np.array(FCSs), np.array(TR_array) def estimate_FCS(self, time_series):