Impaired behavioral inhibition in Fmr1 KO mice is linked to disrupted visual cortex theta oscillations.

Code ↔ Paper

20 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 20 matches

[1] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Local field potential (LFP) analysis ↔ LFP/Neuropixels_LFP/2_reward_LFP-novel.ipynb, lines 162–194 · score 0.94 · 1.05–1.15 s, 0.55–0.7 s, 0.71–0.85 s, 0.86–1 s, 0.86 s, 0.71 s
[2] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Local field potential (LFP) analysis ↔ LFP/V1_HPC/mz_LFP_functions.py, lines 105–230 · score 0.89 · Morlet wavelet, wavelet convolution, phase locking, frequency precision, 3–10, dB
[3] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Classification analyses ↔ ML Code/accvsunits.py, lines 20–119 · score 0.80 · cross validated, balanced accuracy, logistic regression, fold, classifier, bins
[4] § RESULTS › Reduced local field potential dynamics in V1 across conditions ↔ LFP/64Chs_LFP/64Ch_lfp.ipynb, lines 612–655 · score 0.79 · 30–40 Hz, 50–70 Hz, 12–30 Hz, 30–70 Hz, 8–12 Hz, beta
[5] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Classification analyses ↔ ML Code/accvsunits.py, lines 20–119 · score 0.77 · cross validated, balanced accuracy, LBFGS, L2, fold, classification
[6] § RESULTS › Reduced local field potential dynamics in V1 across conditions ↔ LFP/Neuropixels_LFP/2_reward_LFP-novel.ipynb, lines 162–194 · score 0.76 · 0.55–0.7 s, 0.71–0.85 s, 0.86–1 s, 0.86 s, 0.71 s, 0.55 s
[7] § RESULTS › Reduced local field potential dynamics in V1 across conditions ↔ LFP/V1_HPC/plv_phd_afterSaving.ipynb, lines 330–394 · score 0.76 · 30–40 Hz, 50–70 Hz, 12–30 Hz, 30–70 Hz, 8–12 Hz, beta
[8] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Associating theta oscillations with behavior outcomes ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 1577–1636 · score 0.71 · linear regression, standard deviation, theta power, linregress, fit, lick
[9] § RESULTS › Inter-region coupling is decreased in FX while intra-region coherence matches WT ↔ LFP/V1_HPC/plv_phd_afterSaving.ipynb, lines 396–429 · score 0.71 · 0.7–1.5 s, 0.5–0.7 s, 30–70 Hz, V1 HPC, PLV, PHD
[10] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Single unit analysis ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 1294–1336 · score 0.70 · phase domain, LFP phase, Hilbert, SPC, coherence, filter
[11] § RESULTS › Reduced local field potential dynamics in V1 across conditions ↔ LFP/V1_HPC/2_HPC_TF.ipynb, lines 326–358 · score 0.69 · 30–40 Hz, 50–70 Hz, 12–30 Hz, 30–70 Hz, 8–12 Hz, frequency band
[12] § RESULTS › Reduced local field potential dynamics in V1 across conditions ↔ LFP/Neuropixels_LFP/3_TF_analysis-novel.ipynb, lines 241–262 · score 0.67 · 30–40 Hz, 50–70 Hz, 12–30 Hz, 30–70 Hz, 8–12 Hz, 4–8 Hz
[13] § RESULTS › Correlating behavior outcomes with theta oscillation activity ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 1577–1636 · score 0.64 · Linear regressions, Incorrect trials, theta power, fit, correlation, licking
[14] § STAR★METHODS › QUANTIFICATION AND STATISTICAL ANALYSIS ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 1638–1654 · score 0.61 · Spearman correlation coefficient, permutation, licking, behavior, theta, HPC
[15] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Classification analyses ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 781–810 · score 0.60 · logistic regression model, accuracy score, trained, stimulus
[16] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Single unit analysis ↔ LFP/V1_HPC/1_HPC_LFP.ipynb, lines 247–291 · score 0.60 · 30–70 Hz, 12–30, 8–12, 2–4, fft, AUC
[17] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Single unit analysis ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 355–393 · score 0.60 · 30–70 Hz, 12–30, 8–12, 2–4, fft, AUC
[18] § RESULTS › Correlating behavior outcomes with theta oscillation activity ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 1178–1265 · score 0.57 · unit raster, lick raster, LFP traces, CR, Correlating, theta
[19] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Classification analyses ↔ Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, lines 716–779 · score 0.57 · logistic regression model, Go trials, probability, trained
[20] § STAR★METHODS › METHOD DETAILS › Pharmacological interventions › Data acquisition and preprocessing ↔ LFP/64Chs_LFP/64Ch_lfp.ipynb, lines 192–263 · score 0.56 · bandpass filtered, notch, noise, preprocessing, Scipy, electrodes

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 paper is loaded by your browser, with JavaScript: read it at doi.org or on Europe PMC.

The authors' code

Jupyter notebook · 1,672 lines · 79 KB · no license · 8 matches

# %%
import os
import glob
import numpy as np
import random
import pandas as pd
import pickle
import matplotlib as mpl
import matplotlib.pyplot as plt
import seaborn as sns
from matplotlib.patches import Ellipse
import matplotlib.transforms as transforms
import scipy.signal as ssig
from sklearn.metrics import auc
import itertools
import scipy.ndimage as scnd
import scipy.stats as sstat
from statsmodels.formula.api import mixedlm
import mz_LFP_functions as mz_LFP
import mz_ephys_unit_analysis as mz_ena
# %%
mpl.rcParams['pdf.fonttype'] = 42
mpl.rcParams['font.sans-serif']=['Arial', 'Helvetica','Bitstream Vera Sans', 'DejaVu Sans', 'Lucida Grande',
'Verdana', 'Geneva', 'Lucid', 'Avant Garde', 'sans-serif']
rc_pub={'font.size': 20, 'axes.labelsize': 20, 'legend.fontsize': 20,
'axes.titlesize': 25, 'xtick.labelsize': 20, 'ytick.labelsize': 20,
'axes.linewidth':1.5, 'lines.linewidth': 2.0,
'xtick.color': 'black', 'ytick.color': 'black', 'axes.edgecolor': 'black',
'axes.labelcolor':'black','text.color':'black'}
# for publication quality plots
def set_pub_plots(pal=sns.blend_palette(['cyan', 'magenta','gray','crimson','purple'], 5)):
sns.set_style("white")
sns.set_palette(pal)
sns.set_context("poster", font_scale=1.5, rc=rc_pub)
sns.set_style("ticks", {"xtick.major.size": 5, "ytick.major.size": 5})
# to restore the defaults, call plt.rcdefaults()
set_pub_plots()
# %% [markdown]
# # Base variables
# This is primarily base values for ephys data
# %%
insert_depth = 3100 #change this as appropriate
sp_bw_ch = 20/2
surface_ch = np.round(insert_depth/sp_bw_ch)
V1_hip_ch = np.round((insert_depth-1000)/sp_bw_ch)
Hip_thal_ch = np.round((insert_depth-1000-1200)/sp_bw_ch)
CA1_DG_ch = np.round((insert_depth-1000-600)/sp_bw_ch)
samples_tr = 7350 #this is based on the shortest #samples in a trial
sr = 2500
n_chan = 384
rec_length = 3.0 #how long is the arduino triggered
v1_ch = 260
hpc_ch = 110
# %% [markdown]
# # Load ephys data
# output
# - all_trial_lfp <- dataframe of the LFPs for each trial with a 2-d array (#ch x #samp)
# - vr_units <- dataframe of all vr units psth
# - units_fft <- dataframe of all vr units power spectrum analysis
# %% [markdown]
# ## LFP
# %%
all_trial_lfp = pd.read_pickle(r'G:\Neuropixels\02_wtfx_behavior\all_trials.pkl')
display(all_trial_lfp.head())
print(all_trial_lfp.group.unique(), all_trial_lfp.stim_id.unique())
print(all_trial_lfp.groupby('group')['et'].nunique())
# %% [markdown]
# ## Units
# %% [markdown]
# ### PSTH dataframe
# %%
final_df = pd.read_parquet(r"U:\Papers\FX Behavior paper\data\final_V1HPC_OperantNovel_psth.parquet")
vr_units = final_df[final_df.visRes == 'yes']
vr_units.head()
# %%
print(vr_units.region.unique(), vr_units.group.unique(), vr_units.stim_id.unique())
print(vr_units.groupby('group')['et'].nunique())
# %% [markdown]
# ### Create units FFT dataframe
# %%
def get_unit_fft(unit_df):
unit_arr = unit_df[(unit_df.times>0.5)&(unit_df.times<1.8)].zscore.values
freq = np.arange(unit_arr.shape[0]) / unit_arr.shape[0] * 100
freq = freq[:freq.shape[0]//2]
f = np.fft.fft(unit_arr)
magnitude_spectrum = (np.abs(f)[:freq.shape[0]])
return freq, magnitude_spectrum
# %%
units_fft = []
for unit,df in vr_units.groupby(['stim_id','cuid']):
fft_freq, fft_amp = get_unit_fft(df)
bands = [[2,4],[4,8],[8,12],[12,30],[30,70]]
auc_val,band_val = [],[]
for ran in bands:
lower = fft_freq.searchsorted(ran[0], 'left')
upper = fft_freq.searchsorted(ran[1], 'right') -1
val = auc(fft_freq[lower:upper],fft_amp[lower:upper])
auc_val.append(val)
band_val.append(str(ran))
units_fft.append(pd.DataFrame({'band': band_val, 'auc': auc_val, 'cuid': unit[1],
'stim': unit[0], 'group': df.group.unique()[0],
'region': df.region.unique()[0], 'et': df.et.unique()[0]}))
units_fft = pd.concat(units_fft)
# %% [markdown]
# # Load behavior data
# output
# - behavior <- dataframe of trials and behavior label for each mouse
# %%
behavior = pd.read_pickle(r"G:\Neuropixels\02_wtfx_behavior\lick_behavior_rec.pkl")
#replace et CC067431_HP3 with CC067432_HP3 (*this was mislabeled originally)
behavior.loc[behavior["et"] == "CC067431_HP3", "et"] = "CC067432_HP3"
#the stim_id col currently on the df is off, update with true stim_id below
behavior = behavior.drop('stim_id', axis=1)
# %%
# Psuedo random presentation of the stimuli - 25 per row * 6 rows = 150 trials
# 0 -- drifting grating, rewarded stimulus, 100 trials
# 1 -- pink noise, unrewarded stimulus, 50 trials
stim_order = [0,1,1,0,0,1,0,0,1,0,1,0,0,1,0,0,0,1,0,0,0,0,1,0,0,
0,1,0,1,0,1,0,0,1,1,0,0,0,1,0,1,0,0,1,1,0,0,0,0,1,
0,0,0,1,0,0,0,1,0,0,1,0,0,1,0,0,1,1,0,1,0,0,0,0,1,
0,0,0,1,0,1,0,0,1,0,0,1,0,0,0,1,0,0,0,1,0,0,0,0,1,
1,0,0,0,1,0,1,0,0,0,1,0,1,0,0,1,0,0,0,1,0,1,1,0,0,
0,0,1,0,0,1,0,0,1,0,1,0,0,0,1,0,0,1,0,0,1,0,0,1,0]
# Psuedo random distribution of water to the rewarded stimuli
# 0 -- water given -- 80 times
# 1 -- no water given -- 20 times
rew_order = [0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,1,0,1,0,0,0,0,0,1,
0,0,1,0,0,0,0,0,1,0,0,1,0,0,0,0,1,0,0,0,0,1,0,0,0,
0,0,0,1,0,0,0,0,0,1,0,0,0,1,0,0,0,1,0,0,0,1,0,0,0,
0,1,0,0,0,0,0,1,0,0,0,0,0,1,0,0,1,0,0,0,0,1,0,0,0]
overall_order, i = {}, 0
for idx,val in enumerate(stim_order):
if val == 0:
if rew_order[i] == 0:
overall_order[idx]='Go+' #rew
elif rew_order[i] == 1:
overall_order[idx]='Go-' #rew2
i+=1
elif val == 1:
overall_order[idx]='No-Go' #unrew
# %%
#update the stim_id column with the overall order dictionary made above
display(behavior.head(2))
behavior['stim_id'] = behavior['true_tr'].map(overall_order)
display(behavior.head(2))
# %%
#add behavior outcome for each trial based on the stimulus and if licked 2+ times from [0, 1.2] sec
beh_ls = []
for d,dd in behavior.groupby(['et','trial']):
lt = dd.lick_time.values
did_lick = True if sum((lt>0)&(lt<1.2))>=2 else False #checks licks that are between [0, 1.2] time window
stim = dd.stim_id.unique()[0]
if stim == 'Go+':
beh_ls.append(['H']*dd.shape[0] if did_lick else ['M']*dd.shape[0])
elif stim == 'Go-':
beh_ls.append(['H']*dd.shape[0] if did_lick else ['M']*dd.shape[0])
else:
beh_ls.append(['FA']*dd.shape[0] if did_lick else ['CR']*dd.shape[0])
behavior['behav'] = list(itertools.chain(*beh_ls)) #itertools flattens the list of lists to add to the dataframe
print(behavior.behav.unique())
behavior.head()
# %% [markdown]
# # Plot the Results (units and behavior)
# %% [markdown]
# ## Combine ephys and behavior
# %%
mouse_beh = []
for d,dd in behavior.groupby(['et','stim_id','behav']):
counts = dd.true_tr.nunique()
tmp_df=pd.DataFrame({'et':[d[0]],
'stim_id':[d[1]],
'behav':[d[2]],
'counts':counts
})
mouse_beh.append(tmp_df)
mouse_beh = pd.concat(mouse_beh, ignore_index=True)
mouse_beh.head(3)
# %%
et_beh_dict = {}
for d,dd in mouse_beh.groupby('et'):
try:
nmH = dd[dd.behav=='H']['counts'].values[0]
except:
nmH = 0
try:
nmM = dd[dd.behav=='M']['counts'].values[0]
except:
nmM = 0
try:
nmCR = dd[dd.behav=='CR']['counts'].values[0]
except:
nmCR = 0
try:
nmFA = dd[dd.behav=='FA']['counts'].values[0]
except:
nmFA = 0
et_beh_dict[d] = [nmH, nmM, nmCR, nmFA] #build dict like this "{et: [nmH, nmM, nmCR, nmFA]}" to map to below df
# %%
# map dict above with the df above above to give 4 new columns (nmH, nmM, nmCR, nmFA)
units_fft['nmH'] = units_fft.et.map(lambda x: et_beh_dict[x][0])
units_fft['nmM'] = units_fft.et.map(lambda x: et_beh_dict[x][1])
units_fft['nmCR'] = units_fft.et.map(lambda x: et_beh_dict[x][2])
units_fft['nmFA'] = units_fft.et.map(lambda x: et_beh_dict[x][3])
# %%
# take the mean AUC across the units for each mouse
theta_fft = units_fft[(units_fft.band=='[4, 8]')]
theta_fft = theta_fft[(theta_fft.stim==0) | (theta_fft.stim==2)]
mean_theta = []
for d,dd in theta_fft.groupby(['group','et','stim','region','band']):
mean_theta.append({'mean_auc':dd.auc.mean(), 'median_auc':dd.auc.median(),
'group':d[0], 'et':d[1], 'stim':d[2], 'region':d[3], 'band':d[4],
'nmH':dd.nmH.unique()[0], 'nmM':dd.nmM.unique()[0], 'nmCR':dd.nmCR.unique()[0], 'nmFA':dd.nmFA.unique()[0]})
mean_theta = pd.DataFrame(mean_theta)
mean_theta.head(3)
# %% [markdown]
# ## Plot it
# %%
y_choice = 'mean_auc' # mean_auc & median_auc
for x_choice in ["nmH", "nmM", "nmCR", "nmFA"]:
sns.relplot(data=mean_theta, x=x_choice, y=y_choice, hue='group', col='stim', row='region',
palette={'WT':'cyan', 'FX':'magenta'}, height=4, aspect=1.2)
plt.show()
for d,dd in mean_theta.groupby('stim'):
print(f"~~~~~~~~~~ {d} ~~~~~~~~~~")
sns.relplot(data=dd, x='nmH', y='mean_auc',
col='behav', row='region',
hue='group', palette={'WT':'cyan', 'FX':'magenta'},
height=3, facet_kws={'sharey': True, 'sharex': False})
plt.suptitle(d)
plt.show()
# %% [markdown]
# # Plot the results (LFP and behavior)
# %%
display(all_trial_lfp.head(2))
display(behavior.head(2))
# %% [markdown]
# ## LFP traces and licks (all mice all trials)
# %%
v1_ch = 260
hpc_ch = 110
show_plots = False
lfp_behav, i = [], 0
# plot the single trial veps with the licks overlaid
for d,dd in all_trial_lfp.groupby(['et','trial']):
#gets the v1 lfp trace for the et/trial pairing
tr_stim = dd.stim_id.unique()[0]
plt_stim = 'Go+' if tr_stim=='0' else ('Go-' if tr_stim=='1' else 'No-Go')
tr_lfp = dd.lfp_data.values[0]
v1_lfp = tr_lfp[v1_ch,:]
hpc_lfp = tr_lfp[hpc_ch,:]
#trims the behavior df to the et/trial pairing
tmp0 = behavior[(behavior.et==d[0])&(behavior.true_tr==d[1])&(behavior.lick_time>-0.5)&(behavior.lick_time<2.5)]
tmp = tmp0['lick_time'].values+0.5 #lick times now zeroed to start of recording
num_licks = len([x for x in tmp if (x>1.0)&(x<1.7)]) #counting the number of licks within the delay period
tr_beh = tmp0.behav.unique()[0] if tmp0.size else ('CR' if plt_stim=='No-Go' else 'M') #OR else 'nolick'
#save to new df --- columns = trial, et, group, stim_id, v1_lfp, hpc_lfp, num_licks, behav
lfp_behav.append([d[1], d[0], dd.group.unique()[0], plt_stim, v1_lfp.astype('float16'), hpc_lfp.astype('float16'), num_licks, tr_beh])
#plot the lfp traces and the licks overlaid
if show_plots==False:
if i%500 == 0:
print(f"Done with {i} out of {all_trial_lfp.et.nunique()*all_trial_lfp.trial.nunique()}") #loading bar if not plotting LFPs and lick rasters
i+=1
else:
x_times = np.linspace(0,3,v1_lfp.shape[0])
plt_color = 'cyan' if dd.group.unique()[0]=='WT' else 'magenta'
fig,ax = plt.subplots(1,2, sharex=True, figsize=(15,2.5))
ax[0].axvspan(0.5,0.7, color='grey', alpha=0.2)
ax[1].axvspan(0.5,0.7, color='grey', alpha=0.2)
if plt_stim=='Go+':
ax[0].axvline(1.7, color='royalblue')
ax[1].axvline(1.7, color='royalblue')
elif plt_stim=='Go-':
ax[0].axvline(1.7, color='grey')
ax[1].axvline(1.7, color='grey')
else:
ax[0].axvline(1.7, color='grey', ls='--')
ax[1].axvline(1.7, color='grey', ls='--')
ax[0].plot(x_times, scnd.gaussian_filter1d(v1_lfp, sigma=20), color=plt_color) #plot gaussian filtered V1 vep
ax[0].scatter(x=tmp, y=[-500]*len(tmp), marker='|', c='black', s=70) #plot licks
ax[1].plot(x_times, scnd.gaussian_filter1d(hpc_lfp, sigma=20), color=plt_color) #plot gaussian filtered V1 vep
ax[1].scatter(x=tmp, y=[-800]*len(tmp), marker='|', c='black', s=70) #plot licks
ax[0].set_xlabel('Time (s)')
ax[1].set_xlabel('Time (s)')
ax[0].set_ylabel('V1')
ax[0].set_ylim([-600,600])
ax[1].set_ylim([-900,900])
ax[0].set_title(f"{d[0]} ~~ Trial {d[1]} ~~ {plt_stim} ~~ {tr_beh}")
ax[1].set_ylabel('HPC')
sns.despine()
plt.show()
lfp_behav = pd.DataFrame(lfp_behav, columns=['trial','et','group','stim_id','v1_lfp','hpc_lfp','num_licks','behav'])
# %% [markdown]
# ## LFP oscillation FFT and behavior
# %%
def get_lfp_fft(lfp_arr, sr=2500, limit_time=False):
if limit_time:
time_cut_arr = lfp_arr[int(sr*limit_time[0]):int(sr*limit_time[1])] #limits the LFP trace to [0.5, 1.7] and quantify power on that time
else:
time_cut_arr = lfp_arr
freq = np.arange(time_cut_arr.shape[0]) / time_cut_arr.shape[0] * sr #sr of the lfp traces
freq = freq[:freq.shape[0]//2]
f = np.fft.fft(time_cut_arr)
magnitude_spectrum = (np.abs(f)[:freq.shape[0]])
return freq, magnitude_spectrum
# %%
lfp_fft = []
for lfp,df in lfp_behav.groupby(['et','trial']):
# for lfp,df in lfp_behav[lfp_behav.trial<100].groupby(['et','trial']): #limit the # of trials to account for loss of motivation in later ones? ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
fft_freq, fft_amp = get_lfp_fft(df.v1_lfp.values[0], limit_time=[0.5,1.8])
bands = [[2,4],[4,8],[8,12],[12,30],[30,70]]
auc_val,band_val,mean_ls,peak_val = [],[],[],[]
for ran in bands:
lower = fft_freq.searchsorted(ran[0], 'left')
upper = fft_freq.searchsorted(ran[1], 'right') -1
val = auc(fft_freq[lower:upper],fft_amp[lower:upper])
mean_val = np.mean(fft_amp[lower:upper])
peak = np.max(fft_amp[lower:upper])
auc_val.append(val/100000) #scale factor for the auc values
mean_ls.append(mean_val/10000)
peak_val.append(peak/10000)
band_val.append(str(ran))
lfp_fft.append(pd.DataFrame({'trial':lfp[1], 'band':band_val, 'mean_val':mean_ls, 'peak_val':peak_val, 'auc':auc_val, 'et':lfp[0], 'region':'V1',
'group':df.group.unique()[0], 'stim_id':df.stim_id.unique()[0], 'num_licks':df.num_licks.unique()[0], 'behav':df.behav.unique()[0]}))
#repeat above for the HPC channel
fft_freq, fft_amp = get_lfp_fft(df.hpc_lfp.values[0], limit_time=[0.5,1.8])
auc_val,band_val,mean_ls,peak_val = [],[],[],[]
for ran in bands:
lower = fft_freq.searchsorted(ran[0], 'left')
upper = fft_freq.searchsorted(ran[1], 'right') -1
val = auc(fft_freq[lower:upper],fft_amp[lower:upper])
mean_val = np.mean(fft_amp[lower:upper])
peak = np.max(fft_amp[lower:upper])
auc_val.append(val/100000)
mean_ls.append(mean_val/10000)
peak_val.append(peak/10000)
band_val.append(str(ran))
lfp_fft.append(pd.DataFrame({'trial':lfp[1], 'band':band_val, 'mean_val':mean_ls, 'peak_val':peak_val, 'auc':auc_val, 'et':lfp[0], 'region':'HPC',
'group':df.group.unique()[0], 'stim_id':df.stim_id.unique()[0], 'num_licks':df.num_licks.unique()[0], 'behav':df.behav.unique()[0]}))
lfp_fft = pd.concat(lfp_fft)
# replace instances of H2 & M2 with H & M
lfp_fft['behav'] = lfp_fft['behav'].replace(['H2','M2'], ['H','M'])
lfp_fft.head()
# %% [markdown]
# ## Plotting results of all trials
# N = #trials
# %%
do_plot=False
feature_choice = 'peak_val'
for d,dd in lfp_fft.groupby(['stim_id','region']):
if do_plot:
fig,ax = plt.subplots(1,2, sharex=True, figsize=(15,3.5))
plt_order = ['CR','FA'] if d[0]=='No-Go' else ['H','M']
print("2-sided KS test comparing the distributions of WT & FX")
for i,b in enumerate(plt_order):
plt_df = dd[(dd.band=="[4, 8]")&(dd.behav==b)]
x,y = plt_df[plt_df.group=='WT'][feature_choice].values, plt_df[plt_df.group=='FX'][feature_choice].values
res = sstat.ks_2samp(x,y)
plt_stats = "***" if res.pvalue < 0.001 else ("**" if res.pvalue<0.01 else ("*" if res.pvalue <0.05 else "ns"))
if do_plot:
sns.histplot(data=plt_df, x=feature_choice, binwidth=1,
hue='group', hue_order=['WT','FX'],
stat='density', common_norm=False,
kde=True, legend=False, ax=ax[i])
ax[i].set_title(f"{d} ~~ {b}")
if feature_choice == 'peak_val':
ax[i].text(x=25, y=0.06, s=plt_stats)
else:
ax[i].text(x=6, y=0.2, s=plt_stats)
print(f"{d} {b} -- N: WT={len(x)}, F={len(y)} -- p={res.pvalue} --", plt_stats)
if do_plot:
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUChistogram_{d}.pdf", transparent=True)
plt.show()
# %%
for d,dd in lfp_fft.groupby('stim_id'):
aver_df = dd[dd.band=="[4, 8]"]
plt_order = ['CR','FA'] if d=='No-Go' else ['H','M']
sns.catplot(data=aver_df, x='behav', y='mean_val', kind='bar', col='region', order=plt_order,
hue='group', hue_order=['WT','FX'], palette={'WT':'cyan','FX':'magenta'},
height=3.5, aspect=1)
# plt.ylim([0,6.5])
plt.suptitle(d)
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUCbarplot_{d}.pdf", transparent=True)
plt.show()
# %%
feature_choice = 'peak_val'
#Mann Whitney U to compare pairwise within each band for each region
print('~~~~~~~~~~~~ Comparing WT/FX for each region/stim/behav ~~~~~~~~~~~~')
for d,dd in lfp_fft[lfp_fft.band=="[4, 8]"].groupby(['region','stim_id','behav']):
x,y = dd[dd.group=='WT'][feature_choice].values, dd[dd.group=='FX'][feature_choice].values
U,pval = sstat.mannwhitneyu(x,y)
print(f"{d} -- N (#trials, all mice): WT={len(x)}, FX={len(y)} -- pvals: {pval} -- ",
'***' if pval<0.001 else ('**' if pval<0.01 else ('*' if pval<0.05 else 'ns')))
#Mann Whitney U to compare pairwise within each group for each region
print('\n~~~~~~~~~~~~ Comparing behav for each region/stim/group ~~~~~~~~~~~~')
for d,dd in lfp_fft[lfp_fft.band=="[4, 8]"].groupby(['region','stim_id','group']):
beh = dd.behav.unique()
x,y = dd[dd.behav==beh[0]][feature_choice].values, dd[dd.behav==beh[1]][feature_choice].values
U,pval = sstat.mannwhitneyu(x,y)
print(f"{d} -- N (#trials, all mice): WT={len(x)}, FX={len(y)} -- pvals: {pval} -- ",
'***' if pval<0.001 else ('**' if pval<0.01 else ('*' if pval<0.05 else 'ns')))
# %% [markdown]
# # Remake df, grouped by mouse
# This is for each mouse, meaning my N = #mice NOT #trials like above
# %%
trial_lims = [0,150]
lfp_avg_fft = []
for d,dd in lfp_fft[(lfp_fft.trial>=trial_lims[0])&(lfp_fft.trial<=trial_lims[1])].groupby(['et','region','band','group','stim_id','behav']):
lfp_avg_fft.append(pd.DataFrame({'band':d[2], 'auc':dd.auc.mean(), 'peak_val':dd.peak_val.mean(), 'mean_val':dd.mean_val.mean(),
'et':d[0], 'region':d[1], 'group':d[3], 'stim_id':d[4], 'behav':d[5]}, index=[0]))
lfp_avg_fft = pd.concat(lfp_avg_fft, ignore_index=True)
lfp_avg_fft.head()
# %% [markdown]
# # other ideas...
# %%
# remove mice that are inactive (less than n HIT trials)
active_ets = []
for d,dd in lfp_fft[lfp_fft.band=='[4, 8]'].groupby('et'):
num_Hits = dd.groupby('behav').trial.nunique().to_dict()['H'] if 'H' in dd.groupby('behav').trial.nunique().to_dict() else 0
if num_Hits>=5:
active_ets.append(d)
active_ets.remove("CC067489_HP2") #removing an outlier found from the plots below
print(active_ets)
# %% [markdown]
# # Previous plot - histogram
# %%
feature_choice = 'auc' #'auc', 'mean_val', 'peak_val'
for d,dd in lfp_avg_fft[lfp_avg_fft.et.isin(active_ets)].groupby('stim_id'):
plt_data = dd[dd.band=='[4, 8]']
sns.displot(data=plt_data, x=feature_choice, hue='group', hue_order=['WT','FX'],
kde=True, rug=False, stat='density', binwidth=0.35, col='behav', row='region', height=4, aspect=1.2)
plt.ylim([0,0.26])
plt.xlim([0,5])
plt.suptitle(d)
print("2-sided KS test comparing the distributions of WT & FX")
for e,ee in plt_data.groupby(['region','behav']):
x,y = ee[ee.group=='WT'][feature_choice].values, ee[ee.group=='FX'][feature_choice].values # N = #mice bc behavior outcome is separate (ie Go+ has Hit and Miss)
res = sstat.ks_2samp(x,y)
plt_stats = "***" if res.pvalue < 0.001 else ("**" if res.pvalue<0.01 else ("*" if res.pvalue <0.05 else "ns"))
print(f"{d} {e} -- N: WT={len(x)}, F={len(y)} -- p={res.pvalue} --", plt_stats)
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUChistogram_{d}.pdf", transparent=True)
plt.show()
# %%
feature_choice = 'auc' #'auc', 'mean_val', 'peak_val'
for d,dd in lfp_avg_fft[lfp_avg_fft.et.isin(active_ets)].groupby('stim_id'):
plt_data = dd[dd.band=='[4, 8]']
sns.displot(data=plt_data, x=feature_choice, hue='behav', palette={'H':'royalblue','M':'grey','CR':'royalblue','FA':'grey'},
kde=True, rug=False, stat='density', binwidth=0.35, col='group', row='region', height=4, aspect=1.2)
plt.ylim([0,0.26])
plt.xlim([0,5])
plt.suptitle(d)
print("2-sided KS test comparing the distributions of behavior within group")
for e,ee in plt_data.groupby(['region','group']):
if (d=='Go+')|(d=='Go-'):
x,y = ee[ee.behav=='H'][feature_choice].values, ee[ee.behav=='M'][feature_choice].values # N = #mice bc behavior outcome is separate (ie Go+ has Hit and Miss)
else:
x,y = ee[ee.behav=='CR'][feature_choice].values, ee[ee.behav=='FA'][feature_choice].values
res = sstat.ks_2samp(x,y)
plt_stats = "***" if res.pvalue < 0.001 else ("**" if res.pvalue<0.01 else ("*" if res.pvalue <0.05 else "ns"))
print(f"{d} {e} -- N: WT={len(x)}, F={len(y)} -- p={res.pvalue} --", plt_stats)
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUChistogram_{d}.pdf", transparent=True)
plt.show()
# %% [markdown]
# # ! Previous plot - scatter (trial # and AUC, split by behavior)
# %%
feature_choice = 'auc' #'auc', 'mean_val', 'peak_val'
for d,dd in lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et.isin(active_ets))].groupby(['region','stim_id']):
behav_ls = dd.behav.unique()
fig,ax=plt.subplots(1,2,figsize=(14,3), sharey=True)
wt_beh_thresh = dd[dd.group=='WT'][feature_choice].mean() + dd[dd.group=='WT'][feature_choice].std() #mean + 1std to use as the threshold for WT
fx_beh_thresh = dd[dd.group=='FX'][feature_choice].mean() + dd[dd.group=='FX'][feature_choice].std() #mean + 1std to use as the threshold for FX
ax[0].hlines(y=[wt_beh_thresh, fx_beh_thresh], xmin=0, xmax=150, colors=['cyan','magenta'], ls='--')
ax[1].hlines(y=[wt_beh_thresh, fx_beh_thresh], xmin=0, xmax=150, colors=['cyan','magenta'], ls='--')
sns.scatterplot(data=dd[dd.behav==behav_ls[0]], x='trial', y=feature_choice, s=50, hue='group', palette={'WT':'cyan', 'FX':'magenta'}, ax=ax[0], legend=False)
sns.scatterplot(data=dd[dd.behav==behav_ls[1]], x='trial', y=feature_choice, s=50, hue='group', palette={'WT':'cyan', 'FX':'magenta'}, ax=ax[1], legend=False)
ax[0].set_title(f'{d} -- {behav_ls[0]}')
ax[1].set_title(f'{d} -- {behav_ls[1]}')
ax[0].set_xlim([-1,151])
ax[1].set_xlim([-1,151])
ax[0].set_ylim([0,6.5])
ax[1].set_ylim([0,6.5])
ax[0].set_xlabel('Trial')
ax[0].set_ylabel('AUC (a.u.)')
ax[1].set_xlabel('Trial')
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUCscatterplot_{d[0]}_{d[1]}_v2.pdf", transparent=True)
plt.show()
# %%
# make the same threshold for each stimulus group, not each stimulus/behavior group
feature_choice = 'auc' #'auc', 'mean_val', 'peak_val'
group_thresh = True
for d,dd in lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et.isin(active_ets))].groupby(['region','stim_id']):
if d[1] == 'No-Go':
behav_ls = ['CR', 'FA']
else:
behav_ls = ['H', 'M']
fig,ax=plt.subplots(2,2,figsize=(14,6), sharey=True)
ax=ax.flatten()
if group_thresh:
wt_beh_thresh = dd[dd.group=='WT'][feature_choice].mean() + dd[dd.group=='WT'][feature_choice].std() #mean + 1std to use as the threshold for WT
fx_beh_thresh = dd[dd.group=='FX'][feature_choice].mean() + dd[dd.group=='FX'][feature_choice].std() #mean + 1std to use as the threshold for FX
# ax[0].hlines(y=[wt_beh_thresh,fx_beh_thresh], xmin=0, xmax=150, colors=['cyan','magenta'], ls='--')
# ax[1].hlines(y=[wt_beh_thresh,fx_beh_thresh], xmin=0, xmax=150, colors=['cyan','magenta'], ls='--')
ax[0].hlines(y=[wt_beh_thresh], xmin=0, xmax=150, colors=['cyan'], ls='--')
ax[1].hlines(y=[wt_beh_thresh], xmin=0, xmax=150, colors=['cyan'], ls='--')
ax[2].hlines(y=[fx_beh_thresh], xmin=0, xmax=150, colors=['magenta'], ls='--')
ax[3].hlines(y=[fx_beh_thresh], xmin=0, xmax=150, colors=['magenta'], ls='--')
else:
beh_thresh = dd[feature_choice].mean() + (1.5*dd[feature_choice].std()) #mean + 1.5std to use as the threshold for FX
ax[0].axhline(y=beh_thresh, color='grey', ls='--')
ax[1].axhline(y=beh_thresh, color='grey', ls='--')
# sns.scatterplot(data=dd[dd.behav==behav_ls[0]], x='trial', y=feature_choice, s=50, hue='group', palette={'WT':'cyan', 'FX':'magenta'}, ax=ax[0], legend=False)
# sns.scatterplot(data=dd[dd.behav==behav_ls[1]], x='trial', y=feature_choice, s=50, hue='group', palette={'WT':'cyan', 'FX':'magenta'}, ax=ax[1], legend=False)
sns.lineplot(data=dd[(dd.behav==behav_ls[0])&(dd.group=='WT')], x='trial', y=feature_choice, hue='et', palette='tab10', ax=ax[0], legend=False)
sns.lineplot(data=dd[(dd.behav==behav_ls[1])&(dd.group=='WT')], x='trial', y=feature_choice, hue='et', palette='tab10', ax=ax[1], legend=False)
sns.lineplot(data=dd[(dd.behav==behav_ls[0])&(dd.group=='FX')], x='trial', y=feature_choice, hue='et', palette='tab10', ax=ax[2], legend=False)
sns.lineplot(data=dd[(dd.behav==behav_ls[1])&(dd.group=='FX')], x='trial', y=feature_choice, hue='et', palette='tab10', ax=ax[3], legend=False)
ax[0].set_title(f'{d} -- {behav_ls[0]}')
ax[1].set_title(f'{d} -- {behav_ls[1]}')
# ax[0].set_xlim([-1,151])
# ax[1].set_xlim([-1,151])
# ax[0].set_ylim([0,6.5])
# ax[1].set_ylim([0,6.5])
# ax[0].set_xlabel('Trial')
# ax[0].set_ylabel('AUC (a.u.)')
# ax[1].set_xlabel('Trial')
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\allTrial_AUCscatterplot_{d[0]}_{d[1]}.pdf", transparent=True)
plt.show()
# %% [markdown]
# # ! Threshold and find probability of Hit or Miss above mean+1std
# %%
num_std = 1
feature_choice = 'auc' #'auc', 'mean_val', 'peak_val'
above_thresh_df = []
for d,dd in lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et.isin(active_ets))].groupby(['region','stim_id']):
wt_dd, fx_dd = dd[dd.group=='WT'], dd[dd.group=='FX']
wt_beh_thresh = wt_dd[feature_choice].mean() + (num_std*wt_dd[feature_choice].std()) #mean + Xstd to use as the threshold for WT
fx_beh_thresh = fx_dd[feature_choice].mean() + (num_std*fx_dd[feature_choice].std()) #mean + Xstd to use as the threshold for FX
wt_above = wt_dd[wt_dd[feature_choice]>=wt_beh_thresh]
for e,ee in wt_above.groupby('et'):
tot_trls = ee.trial.nunique()
corr_perc = (ee[ee.behav=='CR']['trial'].nunique())/tot_trls*100 if d[1]=='No-Go' else (ee[ee.behav=='H']['trial'].nunique())/tot_trls*100
corr_auc = ee[ee.behav=='CR'][feature_choice].mean() if d[1]=='No-Go' else ee[ee.behav=='H'][feature_choice].mean()
incorr_auc = ee[ee.behav=='FA'][feature_choice].mean() if d[1]=='No-Go' else ee[ee.behav=='M'][feature_choice].mean()
above_thresh_df.append(pd.DataFrame({'region':d[0], 'stim_id':d[1], 'group':ee.group.unique(), 'et':e,
'perc':corr_perc, 'corr_auc':corr_auc, 'incorr_auc':incorr_auc, 'behav':'H/CR'}))
above_thresh_df.append(pd.DataFrame({'region':d[0], 'stim_id':d[1], 'group':ee.group.unique(), 'et':e,
'perc':100-corr_perc, 'corr_auc':corr_auc, 'incorr_auc':incorr_auc, 'behav':'M/FA'}))
fx_above = fx_dd[fx_dd[feature_choice]>=fx_beh_thresh]
for e,ee in fx_above.groupby('et'):
tot_trls = ee.trial.nunique()
corr_perc = (ee[ee.behav=='CR']['trial'].nunique())/tot_trls*100 if d[1]=='No-Go' else (ee[ee.behav=='H']['trial'].nunique())/tot_trls*100
corr_auc = ee[ee.behav=='CR'][feature_choice].mean() if d[1]=='No-Go' else ee[ee.behav=='H'][feature_choice].mean()
incorr_auc = ee[ee.behav=='FA'][feature_choice].mean() if d[1]=='No-Go' else ee[ee.behav=='M'][feature_choice].mean()
above_thresh_df.append(pd.DataFrame({'region':d[0], 'stim_id':d[1], 'group':ee.group.unique(), 'et':e,
'perc':corr_perc, 'corr_auc':corr_auc, 'incorr_auc':incorr_auc, 'behav':'H/CR'}))
above_thresh_df.append(pd.DataFrame({'region':d[0], 'stim_id':d[1], 'group':ee.group.unique(), 'et':e,
'perc':100-corr_perc, 'corr_auc':corr_auc, 'incorr_auc':incorr_auc, 'behav':'M/FA'}))
above_thresh_df = pd.concat(above_thresh_df, ignore_index=True)
above_thresh_df.head()
# %%
num_std = 1.5 # standard deviation you want to have contained within the ellipse
dot_size = 50 # aesthetic choice for the size of the scatter plot dots
plot_linReg = True
pt1_above_thresh_df = above_thresh_df[above_thresh_df.behav=="H/CR"]
idx=0
fig,ax = plt.subplots(1,2, figsize=(10,3), sharey=True, sharex=True)
for d,dd in pt1_above_thresh_df.groupby("region"): # THIS IS INCORPORATING BOTH THE GO AND NO_GO TRIALS TOGETHER FOR THE CORRECT RESPONSES ONLY
dd = dd.dropna()
x = dd[(dd.group=='WT')]['perc'].values
y = dd[(dd.group=='WT')]['corr_auc'].values
x1 = dd[(dd.group=='FX')]['perc'].values
y1 = dd[(dd.group=='FX')]['corr_auc'].values
# plot linear regression of the points & print R**2 value of the fit
if plot_linReg:
res, res1 = sstat.linregress(x, y), sstat.linregress(x1, y1)
print(f"------------- {d} -------------")
print(f'Slope: WT {res.slope:.6f}, FX {res1.slope:.6f}')
print(f'R squared: WT {res.rvalue**2:.6f}, FX {res1.rvalue**2:.6f}')
print(f'pval: WT {res.pvalue:.6f}, FX {res1.pvalue:.6f}')
ax[idx].plot(x, res.intercept+res.slope*x, color='cyan')
ax[idx].plot(x1, res1.intercept+res1.slope*x1, color='magenta')
# Plot scatter
ax[idx].scatter(x, y, s=dot_size, edgecolors=None, color='cyan', label='WT')
ax[idx].scatter(x1, y1, s=dot_size, edgecolors=None, color='magenta', label='FX')
ax[idx].set_xlabel("Correct %")
ax[idx].set_ylabel('Theta AUC')
ax[idx].set_title(d)
idx+=1
sns.despine()
plt.show()
# %%
# the results above show a significant correlation with correct behavior response and Theta AUC for WT V1!! this is in theresholded trials above mean+1std
# %%
g = sns.catplot(data=above_thresh_df, x='behav', y='perc', kind='bar', hue='group', row='region', row_order=['V1', 'HPC'], col='stim_id', errorbar='se', height=3.5, aspect=1.1)
g.set(xlabel='', ylabel='% of Trials')
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\barplot_thresholdTrialPercent.pdf", transparent=True)
plt.show()
print('----------- WT/FX -----------') # stats for comparing WT/FX
for d,dd in above_thresh_df.groupby(["region","stim_id"]):
for e,ee in dd.groupby('behav'):
x,y = ee[ee.group=='WT'].perc.values, ee[ee.group=='FX'].perc.values
num_x,num_y = ee[ee.group=='WT'].et.nunique(), ee[ee.group=='FX'].et.nunique()
u,p = sstat.mannwhitneyu(x,y)
print(d,e,num_x,num_y,p,"***" if p<0.001 else ("**" if p<0.01 else "*" if p<0.05 else "ns"))
print('----------- C/IC -----------') # stats for comparing within group
for d,dd in above_thresh_df.groupby(["region","stim_id"]):
for e,ee in dd.groupby('group'):
x,y = ee[ee.behav=='H/CR'].perc.values, ee[ee.behav=='M/FA'].perc.values
num_x,num_y = ee[ee.behav=='H/CR'].et.nunique(), ee[ee.behav=='M/FA'].et.nunique()
u,p = sstat.mannwhitneyu(x,y)
print(d,e,num_x,num_y,p,"***" if p<0.001 else ("**" if p<0.01 else "*" if p<0.05 else "ns"))
# %% [markdown]
# # Logistic regression on the trial by trial data - MZ 03.07.25
# %%
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score, classification_report, confusion_matrix, roc_curve, auc
import statsmodels.api as sm
theta_lfp_fft = lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et.isin(active_ets))]
theta_lfp_fft.head(2)
# %%
is_plot = True
run_stats = True
for d,dd in theta_lfp_fft.groupby(['group', 'region']): #combine all Go and No-Go trials together
# for d,dd in theta_lfp_fft[theta_lfp_fft.stim_id!='No-Go'].groupby(['group','region']): #Only look at Go trials
# for d,dd in theta_lfp_fft[theta_lfp_fft.stim_id=='No-Go'].groupby(['group','region']): #Only look at No-Go trials
X = dd.auc.to_numpy().reshape(-1,1)
foo = dd.behav.values
y_binary = np.select([foo=='H', foo=='M', foo=='CR', foo=='FA'], [1,0,1,0], foo).astype('int32') #replace behavior with 1 or 0 if correct or incorrect (BOTH GO AND NO-GO INCLUDED)
# Split the data into training and testing sets
X_train, X_test, y_train, y_test = train_test_split(X, y_binary, test_size=0.2, random_state=42)
# Standardize features
scaler = StandardScaler()
X_train = scaler.fit_transform(X_train)
X_test = scaler.transform(X_test)
# Train the Logistic Regression model
model = LogisticRegression()
model.fit(X_train, y_train)
# Evaluate the model
y_pred = model.predict(X_test)
accuracy = accuracy_score(y_test, y_pred)
if run_stats:
X = dd.auc.to_numpy().reshape(-1,1)
foo = dd.behav.values
y_binary = np.select([foo=='H', foo=='M', foo=='CR', foo=='FA'], [1,0,1,0], foo).astype('int32') #replace behavior outcome with 1 or 0
# Split the data into training and testing sets
X_train, X_test, y_train, y_test = train_test_split(X, y_binary, test_size=0.2, random_state=42)
# building the model and fitting the data
log_reg = sm.Logit(y_train, X_train).fit()
pstar = '***' if log_reg.pvalues < 0.001 else ('**' if log_reg.pvalues < 0.01 else ('*' if log_reg.pvalues < 0.05 else 'ns'))
print(log_reg.summary())
print(f"pvals: {log_reg.pvalues} -- {pstar}\n")
if is_plot:
# Plot logistic regression fit Curve
fig,ax = plt.subplots(1,3, figsize=(13,4))
data = pd.DataFrame({'auc': dd.auc.to_numpy(), 'behav': y_binary})
sns.regplot(x=data['auc'], y=data['behav'], data=data, logistic=True, ci=None, color="cyan" if d[0]=='WT' else 'magenta',
scatter_kws={'s': 10}, line_kws={'color': 'gray'}, ax=ax[0])
ax[0].set_title(d)
ax[0].set_xlim([0,9])
sns.despine(ax=ax[0])
# Plot ROC Curve
y_prob = model.predict_proba(X_test)[:, 1]
fpr, tpr, thresholds = roc_curve(y_test, y_prob)
roc_auc = auc(fpr, tpr)
ax[1].plot(fpr, tpr, color='darkorange', lw=2, label=f'ROC Curve (AUC = {roc_auc:.2f})')
ax[1].plot([0, 1], [0, 1], color='grey', lw=2, linestyle='--', label='Random')
ax[1].set_xlabel('False Positive Rate')
ax[1].set_ylabel('True Positive Rate')
ax[1].set_title(f'Accuracy: {accuracy * 100:.2f}%')
sns.despine(ax=ax[1])
# Confusion matrix
cm = confusion_matrix(y_test, y_pred)
print(cm)
sns.heatmap(cm, annot=True, fmt='d', cmap='Blues' if d[0]=='WT' else 'Reds', ax=ax[2])
ax[2].set_xlabel('Predicted')
ax[2].set_ylabel('Actual')
plt.tight_layout()
plt.show()
# %%
# similar code to the above cell, just doesn't plot and saves accuracy values to a DataFrame
accuracy_ls = []
for stim in ['Go', 'No-Go']:
if stim == 'Go':
foo_df = theta_lfp_fft[theta_lfp_fft.stim_id!='No-Go']
else:
foo_df = theta_lfp_fft[theta_lfp_fft.stim_id=='No-Go']
for d,dd in foo_df.groupby(['group', 'region', 'et']): #combine all Go and No-Go trials together
try:
X = dd.auc.to_numpy().reshape(-1,1)
foo = dd.behav.values
y_binary = np.select([foo=='H', foo=='M', foo=='CR', foo=='FA'], [1,0,1,0], foo).astype('int32') #replace behavior with 1 or 0 if correct or incorrect (BOTH GO AND NO-GO INCLUDED)
# Split the data into training and testing sets
X_train, X_test, y_train, y_test = train_test_split(X, y_binary, test_size=0.2, random_state=42)
# Standardize features
scaler = StandardScaler()
X_train = scaler.fit_transform(X_train)
X_test = scaler.transform(X_test)
# Train the Logistic Regression model
model = LogisticRegression()
model.fit(X_train, y_train)
# Evaluate the model
y_pred = model.predict(X_test)
accuracy = accuracy_score(y_test, y_pred)
accuracy_ls.append(pd.DataFrame({'Accuracy': [accuracy], 'et': d[2], 'group': d[0], 'region':d[1], 'stim_id':[stim]}))
except:
continue
accuracy_df = pd.concat(accuracy_ls, ignore_index=True)
accuracy_df.head()
# %%
sns.catplot(data=accuracy_df, kind="bar", x="group", y="Accuracy", order=['WT', 'FX'],
col='stim_id', col_order=['Go', 'No-Go'], row='region', row_order=['V1', 'HPC'],
height=4, aspect=1.25)
for d,dd in accuracy_df.groupby(['stim_id', 'region']):
x = dd[dd.group=='WT'].Accuracy.values
y = dd[dd.group=='FX'].Accuracy.values
res = sstat.shapiro(x)
print("wt --", res)
res = sstat.shapiro(y)
print("fx --", res)
res = sstat.ttest_ind(x, y)
pstar = '***' if res.pvalue < 0.001 else ('**' if res.pvalue < 0.01 else ('*' if res.pvalue < 0.05 else "ns"))
print(f"{d} -- {res} -- {pstar}")
for d,dd in accuracy_df.groupby(['group', 'region']):
x = dd[dd.stim_id=='Go'].Accuracy.values
y = dd[dd.stim_id=='No-Go'].Accuracy.values
res = sstat.ttest_ind(x, y)
pstar = '***' if res.pvalue < 0.001 else ('**' if res.pvalue < 0.01 else ('*' if res.pvalue < 0.05 else "ns"))
print(f"{d} -- {res} -- {pstar}")
plt.show()
# %%
# %%
# %%
# %%
# %%
# %% [markdown]
# # ! 2d scatter plots for lick prob and behavior for the LFP data
# %%
osc_characteristic = 'auc' #'auc', 'mean_val', 'peak_val'
go_evidence = []
for d,dd in lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et.isin(active_ets))].groupby(['et','region']):
go_c = dd[(dd.stim_id=='Go+')&(dd.behav=='H')|(dd.stim_id=='Go-')&(dd.behav=='H')][osc_characteristic].mean()
go_i = dd[(dd.stim_id=='Go+')&(dd.behav=='M')|(dd.stim_id=='Go-')&(dd.behav=='M')][osc_characteristic].mean()
go_all = dd[(dd.stim_id=='Go+')|(dd.stim_id=='Go-')][osc_characteristic].mean()
ng_c = dd[(dd.stim_id=='No-Go')&(dd.behav=='CR')][osc_characteristic].mean()
ng_i = dd[(dd.stim_id=='No-Go')&(dd.behav=='FA')][osc_characteristic].mean()
ng_all = dd[dd.stim_id=='No-Go'][osc_characteristic].mean()
go = (go_c-go_i)/(go_c+go_i)
ng = (ng_c-ng_i)/(ng_c+ng_i)
g_ng = (go-ng)/(go+ng)
go_base = (go_c-go_all)/(dd[(dd.stim_id=='Go+')|(dd.stim_id=='Go-')][osc_characteristic].std())
ng_base = (ng_c-ng_all)/(dd[dd.stim_id=='No-Go'][osc_characteristic].std())
lick_dict = dd.groupby('behav').trial.nunique().to_dict()
d_prime = sstat.norm.ppf(lick_dict['H']/100) - sstat.norm.ppf(lick_dict['FA']/50)
go_evidence.append(pd.DataFrame({'goC':go_c, 'goIC':go_i,'ngC':ng_c, 'ngIC':ng_i, 'go_avg':go, 'ng_avg':ng, 'gng_avg':g_ng, 'go_base':go_base, 'ng_base':ng_base,
'go_licks':lick_dict['H']/100, 'ng_licks':lick_dict['CR']/50, 'gng_licks':(lick_dict['H']+lick_dict['CR'])/150, 'Dprime':d_prime,
'go2_licks':lick_dict['M']/100, 'ng2_licks':lick_dict['FA']/50,
'et':d[0], 'group':dd.group.unique()[0], 'region':d[1]}, index=[0]))
go_evidence = pd.concat(go_evidence, ignore_index=True)
go_evidence.head()
# %%
def confidence_ellipse(x, y, ax, n_std=3.0, facecolor='none', **kwargs):
"""
Create a plot of the covariance confidence ellipse of *x* and *y*.
Parameters
----------
x, y : array-like, shape (n, )
Input data.
ax : matplotlib.axes.Axes
The Axes object to draw the ellipse into.
n_std : float
The number of standard deviations to determine the ellipse's radiuses.
**kwargs
Forwarded to `~matplotlib.patches.Ellipse`
Returns
-------
matplotlib.patches.Ellipse
CODE FROM: https://matplotlib.org/stable/gallery/statistics/confidence_ellipse.html#sphx-glr-gallery-statistics-confidence-ellipse-py
"""
if x.size != y.size:
raise ValueError("x and y must be the same size")
cov = np.cov(x, y)
pearson = cov[0, 1]/np.sqrt(cov[0, 0] * cov[1, 1])
# Using a special case to obtain the eigenvalues of this two-dimensional dataset.
ell_radius_x = np.sqrt(1 + pearson)
ell_radius_y = np.sqrt(1 - pearson)
ellipse = Ellipse((0, 0), width=ell_radius_x * 2, height=ell_radius_y * 2, facecolor=facecolor, **kwargs)
# Calculating the standard deviation of x from the squareroot of the variance and multiplying with the given number of standard deviations.
scale_x = np.sqrt(cov[0, 0]) * n_std
mean_x = np.mean(x)
# calculating the standard deviation of y ...
scale_y = np.sqrt(cov[1, 1]) * n_std
mean_y = np.mean(y)
transf = transforms.Affine2D().rotate_deg(45).scale(scale_x, scale_y).translate(mean_x, mean_y)
ellipse.set_transform(transf + ax.transData)
return ax.add_patch(ellipse)
# %%
num_std = 1.5 # standard deviation you want to have contained within the ellipse
dot_size = 50 # aesthetic choice for the size of the scatter plot dots
plot_linReg = True # adds a linear regression to the 2d scatter (with the ellipse)
# for val in [['go_licks', 'go_base'], ['ng_licks', 'ng_base']]:
# for val in [['go_licks', 'go_avg'], ['ng_licks', 'ng_avg'], ['gng_licks', 'gng_avg'], ['Dprime', 'gng_avg']]:
# for val in [['go_licks', 'goC'], ['ng_licks', 'ngC']]:
for val in [['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]:
my_x, my_y = val[0], val[1]
fig, ax = plt.subplots(1, 2, figsize=(10, 4), sharex=True, sharey=True)
for idx, reg in enumerate(['V1','HPC']):
x = go_evidence[(go_evidence.group=='WT')&(go_evidence.region==reg)][my_x].values
y = go_evidence[(go_evidence.group=='WT')&(go_evidence.region==reg)][my_y].values
x1 = go_evidence[(go_evidence.group=='FX')&(go_evidence.region==reg)][my_x].values
y1 = go_evidence[(go_evidence.group=='FX')&(go_evidence.region==reg)][my_y].values
# plot linear regression of the points & print R**2 value of the fit
if plot_linReg:
res, res1 = sstat.linregress(x, y), sstat.linregress(x1, y1)
print(f"------------- {my_x} -------------")
print(f'{reg} - Slope: WT {res.slope:.6f}, FX {res1.slope:.6f}')
print(f'{reg} - R squared: WT {res.rvalue**2:.6f}, FX {res1.rvalue**2:.6f}')
print(f'{reg} - pval: WT {res.pvalue:.6f}, FX {res1.pvalue:.6f}')
ax[idx].plot(x, res.intercept+res.slope*x, color='cyan')
ax[idx].plot(x1, res1.intercept+res1.slope*x1, color='magenta')
# Plot scatter
ax[idx].scatter(x, y, s=dot_size, edgecolors=None, color='cyan', label='WT')
ax[idx].scatter(x1, y1, s=dot_size, edgecolors=None, color='magenta', label='FX')
# Plot ellipse using above function
confidence_ellipse(x, y, ax[idx], n_std=num_std, edgecolor=None, facecolor='cyan', alpha=0.3)
confidence_ellipse(x1, y1, ax[idx], n_std=num_std, edgecolor=None, facecolor='magenta', alpha=0.3)
if my_x == 'gng_licks':
ax[idx].set_xlim([0,1.1])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("Correct %")
ax[idx].set_ylabel("GNG Theta")
elif my_x == 'Dprime':
ax[idx].set_xlim([-1,2])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("d'")
ax[idx].set_ylabel("GNG Theta")
else:
ax[idx].set_xlim([-0.1,1.1])
ax[idx].set_xlabel("Correct rate")
if (my_y=='go_avg')|(my_y=='ng_avg'):
ax[idx].set_ylabel("Go Theta") if my_x=='go_licks' else ax[idx].set_ylabel("No-Go Theta")
else:
if (my_y=='goIC')|(my_y=='ngIC'):
ax[idx].set_ylabel("Incorrect Go") if my_x=='go2_licks' else ax[idx].set_ylabel("Incorrect No-Go")
ax[idx].set_ylim([0.5,4.5])
ax[idx].set_xlabel("Incorrect rate")
else:
ax[idx].set_ylabel("Correct Go") if my_x=='go_licks' else ax[idx].set_ylabel("Correct No-Go")
ax[idx].set_ylim([0.5,4.5])
ax[idx].set_title(reg)
plt.tight_layout()
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\2dscatter_LFP_incorrectTrials_{my_x}_v1hpc.pdf", transparent=True)
plt.show()
# %%
# Calculate a Spearman correlation coefficient with associated p-value on the 2d scatter plots
def statistic(x): # permute only `x`
return sstat.spearmanr(x, y).statistic
# https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.spearmanr.html <--- Spearman test info
# https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.permutation_test.html <--- Permutation info
# for val in [['go_licks', 'go_avg'], ['ng_licks', 'ng_avg'], ['gng_licks', 'gng_avg'], ['Dprime', 'gng_avg']]:
for val in [['go_licks', 'goC'], ['ng_licks', 'ngC'], ['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]:
my_x, my_y = val[0], val[1]
for idx, reg in enumerate(['V1','HPC']):
for g in ['WT', 'FX']:
x = go_evidence[(go_evidence.group==g)&(go_evidence.region==reg)][my_x].values
y = go_evidence[(go_evidence.group==g)&(go_evidence.region==reg)][my_y].values
res_exact = sstat.permutation_test((x,), statistic, permutation_type='pairings')
print(val, reg, g)
pval = "***" if res_exact.pvalue < 0.001 else ("**" if res_exact.pvalue<0.01 else ("*" if res_exact.pvalue <0.05 else "ns"))
print(f"stat = {res_exact.statistic:.6f} -- pval = {res_exact.pvalue:.6f} -- {pval}")
# %%
# trying a 2-dimensional, 2-sample KS test for comparing the WT and FX groups (https://github.com/syrte/ndtest)
import ndtest
for val in [['go_licks', 'goC'], ['ng_licks', 'ngC'], ['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]: #correct/incorrect trial responses only
my_x, my_y = val[0], val[1]
for idx, reg in enumerate(['V1','HPC']):
x = go_evidence[(go_evidence.group=='WT')&(go_evidence.region==reg)][my_x].values
y = go_evidence[(go_evidence.group=='WT')&(go_evidence.region==reg)][my_y].values
x1 = go_evidence[(go_evidence.group=='FX')&(go_evidence.region==reg)][my_x].values
y1 = go_evidence[(go_evidence.group=='FX')&(go_evidence.region==reg)][my_y].values
P, D = ndtest.ks2d2s(x, y, x1, y1, extra=True)
pval = "***" if P<0.001 else ("**" if P<0.01 else ("*" if P<0.05 else "ns"))
print(f"{my_y} -- {reg} -- d={D:.7f} -- p={P:.7f} -- {pval}")
# %%
# %%
# %%
# %%
# %%
# %%
# # trying a multivariate ANOVA test (https://www.reneshbedre.com/blog/manova-python.html)
# # https://www.statsmodels.org/stable/generated/statsmodels.multivariate.manova.MANOVA.html#statsmodels.multivariate.manova.MANOVA
# from statsmodels.multivariate.manova import MANOVA
# for reg in ['V1','HPC']:
# region_df = go_evidence[go_evidence.region==reg]
# fit = MANOVA.from_formula('go_licks + goC ~ group', data=region_df)
# print(f'goC -- {reg}')
# print(fit.mv_test())
# for reg in ['V1','HPC']:
# region_df = go_evidence[go_evidence.region==reg]
# fit = MANOVA.from_formula('ng_licks + ngC ~ group', data=region_df)
# print(f'goC -- {reg}')
# print(fit.mv_test())
# %%
# # Trying Linear Mixed Model (https://www.statsmodels.org/stable/mixed_linear.html)
# import statsmodels.formula.api as smf
# for reg in ['V1','HPC']:
# region_df = go_evidence[go_evidence.region==reg]
# md = smf.mixedlm("goC ~ go_licks", region_df, groups=region_df["group"], re_formula="~go_licks")
# mdf = md.fit(method=["lbfgs"])
# print(mdf.summary())
# for reg in ['V1','HPC']:
# region_df = go_evidence[go_evidence.region==reg]
# md = smf.mixedlm("ngC ~ ng_licks", region_df, groups=region_df["group"], re_formula="~ng_licks")
# mdf = md.fit(method=["lbfgs"])
# print(mdf.summary())
# %%
# %%
# %% [markdown]
# ## ! 2d scatter plots for lick prob and behavior for the units data
# %%
display(units_fft.head(2))
display(vr_units.head(2))
# %%
osc_characteristic = 'auc'
go_evid_units = []
for d,dd in units_fft[(units_fft.band=="[4, 8]")&(units_fft.et.isin(active_ets))].groupby(['et','region']):
go = dd[(dd.stim==0)|(dd.stim==1)][osc_characteristic].mean()
ng = dd[(dd.stim==2)][osc_characteristic].mean()
g_ng = go/(go+ng)
go_evid_units.append(pd.DataFrame({'go_avg':go, 'ng_avg':ng, 'gng_avg':g_ng, 'et':d[0], 'group':dd.group.unique()[0], 'region':d[1]}, index=[0]))
go_evid_units = pd.concat(go_evid_units, ignore_index=True)
#map lick probability values from LFP df to units df
map_dict = dict(zip(go_evidence['et'], go_evidence['gng_licks']))
go_evid_units['gng_licks'] = go_evid_units['et'].map(map_dict)
map_dict = dict(zip(go_evidence['et'], go_evidence['go_licks']))
go_evid_units['go_licks'] = go_evid_units['et'].map(map_dict)
map_dict = dict(zip(go_evidence['et'], go_evidence['ng_licks']))
go_evid_units['ng_licks'] = go_evid_units['et'].map(map_dict)
map_dict = dict(zip(go_evidence['et'], go_evidence['Dprime']))
go_evid_units['Dprime'] = go_evid_units['et'].map(map_dict)
go_evid_units.head()
# %%
num_std = 1.5 # standard deviation you want to have contained within the ellipse
dot_size = 50 # aesthetic choice for the size of the scatter plot dots
plot_linReg = False # adds a linear regression to the 2d scatter (with the ellipse)
for val in [['go_licks', 'go_avg'], ['ng_licks', 'ng_avg'], ['gng_licks', 'gng_avg'], ['Dprime', 'gng_avg']]:
my_x, my_y = val[0], val[1]
fig, ax = plt.subplots(1, 2, figsize=(10, 4), sharex=True, sharey=True)
for idx, reg in enumerate(['v1','hippo']):
x = go_evid_units[(go_evid_units.group=='WT')&(go_evid_units.region==reg)][my_x].values
y = go_evid_units[(go_evid_units.group=='WT')&(go_evid_units.region==reg)][my_y].values
x1 = go_evid_units[(go_evid_units.group=='FX')&(go_evid_units.region==reg)][my_x].values
y1 = go_evid_units[(go_evid_units.group=='FX')&(go_evid_units.region==reg)][my_y].values
# plot linear regression of the points & print R**2 value of the fit
if plot_linReg:
res, res1 = sstat.linregress(x, y), sstat.linregress(x1, y1)
print(f'{reg} - Slope: WT {res.slope:.6f}, FX {res1.slope:.6f}')
print(f'{reg} - R squared: WT {res.rvalue**2:.6f}, FX {res1.rvalue**2:.6f}')
ax[idx].plot(x, res.intercept+res.slope*x, color='cyan')
ax[idx].plot(x1, res1.intercept+res1.slope*x1, color='magenta')
# Plot scatter
ax[idx].scatter(x, y, s=dot_size, edgecolors=None, color='cyan', label='WT')
ax[idx].scatter(x1, y1, s=dot_size, edgecolors=None, color='magenta', label='FX')
# Plot ellipse using above function
confidence_ellipse(x, y, ax[idx], n_std=num_std, edgecolor='cyan', facecolor='cyan', alpha=0.2)
confidence_ellipse(x1, y1, ax[idx], n_std=num_std, edgecolor='magenta', facecolor='magenta', alpha=0.2)
if my_x == 'gng_licks':
ax[idx].set_xlim([0,1.1])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("Lick probability")
ax[idx].set_ylabel("GNG Theta")
elif my_x == 'Dprime':
ax[idx].set_xlim([-1,2])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("d'")
ax[idx].set_ylabel("GNG Theta")
else:
ax[idx].set_xlim([-0.1,1.1])
# ax[idx].set_ylim([0.4,0.6])
ax[idx].set_xlabel("Lick probability")
ax[idx].set_ylabel("Go Theta") if my_x=='go_licks' else ax[idx].set_ylabel("No-Go Theta")
ax[idx].set_title(f"{reg} units")
plt.tight_layout()
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\2dscatter_units_{my_x}_v1hpc.pdf", transparent=True)
plt.show()
# %% [markdown]
# # ! LFP trace, unit rasters, and lick behavior
# %%
unit_spikes = pd.read_pickle(r"U:\Papers\FX Behavior paper\data\operant_spikes_df.pkl")
unit_spikes.head()
# %%
# #plot rasters for all units to decide which ones to use for the final combo plots in cell below
# my_et = "CC082260_HP2" # WT mouse
# # my_et = "CC082255_HP0" # FX mouse
# stim_id = 0
# plt_color = 'cyan' if my_et == "CC082260_HP2" else 'magenta'
# et_behav = dict(zip(range(150),lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et==my_et)&(lfp_fft.region=='V1')].behav.values))
# et_spk = unit_spikes[unit_spikes.et==my_et]
# et_spk['behav'] = et_spk['trial'].map(et_behav)
# for unit in vr_units[vr_units.et==my_et].cuid.unique():
# print(unit)
# unit_spk = et_spk[et_spk.cuid==unit]
# fig, ax = plt.subplots(2,2, figsize=(15,5), sharex=True, sharey=True)
# ax = ax.flatten()
# ax_titles = ['Hit', "Miss", 'CR', "FA"]
# for idx,axis in enumerate(ax):
# axis.axvspan(0.5 ,0.7, color='grey', alpha=0.2)
# axis.set(ylabel='Trial #')
# axis.set_title(ax_titles[idx])
# #raster plots
# ax[0].plot(unit_spk[unit_spk.behav=='H'].trial_spikes, unit_spk[unit_spk.behav=='H'].trial, '.', ms=8, color=plt_color)
# ax[1].plot(unit_spk[unit_spk.behav=='M'].trial_spikes, unit_spk[unit_spk.behav=='M'].trial, '.', ms=8, color=plt_color)
# ax[2].plot(unit_spk[unit_spk.behav=='CR'].trial_spikes, unit_spk[unit_spk.behav=='CR'].trial, '.', ms=8, color=plt_color)
# ax[3].plot(unit_spk[unit_spk.behav=='FA'].trial_spikes, unit_spk[unit_spk.behav=='FA'].trial, '.', ms=8, color=plt_color)
# sns.despine()
# plt.show()
# %%
# this is dict with {et: [trial #s]}
d = {'CC082260_HP2':[147,56,124,107,122,30], 'CC082255_HP0':[84,129,51,9,119,78]}
wt_v1_units = ['CC082260_HP2_477', 'CC082260_HP2_544', 'CC082260_HP2_567', 'CC082260_HP2_521', 'CC082260_HP2_517',
'CC082260_HP2_548', 'CC082260_HP2_562', 'CC082260_HP2_559']
wt_hpc_units = ['CC082260_HP2_187', 'CC082260_HP2_188', 'CC082260_HP2_201', 'CC082260_HP2_224', 'CC082260_HP2_236',
'CC082260_HP2_239', 'CC082260_HP2_240', 'CC082260_HP2_254']
fx_v1_units = ['CC082255_HP0_564', 'CC082255_HP0_571', 'CC082255_HP0_589', 'CC082255_HP0_577', 'CC082255_HP0_517',
'CC082255_HP0_601', 'CC082255_HP0_554']
fx_hpc_units = ['CC082255_HP0_193', 'CC082255_HP0_212', 'CC082255_HP0_214', 'CC082255_HP0_227', 'CC082255_HP0_259',
'CC082255_HP0_263']
unit_gap = 120
rast_start = 1000
is_saving = False
for et,trs in d.items():
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
et_behav = dict(zip(range(150),lfp_fft[(lfp_fft.band=="[4, 8]")&(lfp_fft.et==et)&(lfp_fft.region=='V1')].behav.values))
et_spk = unit_spikes[unit_spikes.et==et]
et_spk['behav'] = et_spk['trial'].map(et_behav)
if et == "CC082260_HP2":
unit_spk = et_spk[(et_spk.cuid.isin(wt_v1_units))|(et_spk.cuid.isin(wt_hpc_units))] #wt units list for the raster plots
v1_units, hpc_units = wt_v1_units, wt_hpc_units
else:
unit_spk = et_spk[(et_spk.cuid.isin(fx_v1_units))|(et_spk.cuid.isin(fx_hpc_units))]
v1_units, hpc_units = fx_v1_units, fx_hpc_units
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# for plt_tr in trs:
for plt_tr in [7, 20, 146, 101]: #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Make this a list of 0-149 to look for another VEP Mis for Go+ and Go-
dd = all_trial_lfp[(all_trial_lfp.et==et)&(all_trial_lfp.trial==plt_tr)]
#gets the v1 lfp trace for the et/trial pairing
tr_stim = dd.stim_id.unique()[0]
plt_stim = 'Go+' if tr_stim=='0' else ('Go-' if tr_stim=='1' else 'No-Go')
tr_lfp = dd.lfp_data.values[0]
v1_lfp = tr_lfp[v1_ch,:]
hpc_lfp = tr_lfp[hpc_ch,:]
#trims the behavior df to the et/trial pairing
tmp0 = behavior[(behavior.et==et)&(behavior.true_tr==plt_tr)&(behavior.lick_time>-0.5)&(behavior.lick_time<2.5)]
tmp = tmp0['lick_time'].values+0.5 #lick times now zeroed to start of recording
num_licks = len([x for x in tmp if (x>1.0)&(x<1.7)]) #counting the number of licks within the delay period
tr_beh = tmp0.behav.unique()[0] if tmp0.size else ('CR' if plt_stim=='No-Go' else 'M') #OR else 'nolick'
#plot the lfp traces and the licks overlaid
x_times = np.linspace(0,3,v1_lfp.shape[0])
plt_color = 'cyan' if dd.group.unique()[0]=='WT' else 'magenta'
fig,ax = plt.subplots(1,2, sharex=True, sharey=True, figsize=(15,5))
ax[0].axvspan(0.5,0.7, color='grey', alpha=0.2)
ax[1].axvspan(0.5,0.7, color='grey', alpha=0.2)
if plt_stim=='Go+':
ax[0].axvline(1.7, color='royalblue')
ax[1].axvline(1.7, color='royalblue')
elif plt_stim=='Go-':
ax[0].axvline(1.7, color='grey')
ax[1].axvline(1.7, color='grey')
ax[0].plot(x_times, scnd.gaussian_filter1d(v1_lfp, sigma=20), color=plt_color, linewidth=4) #plot gaussian filtered V1 vep
ax[0].scatter(x=tmp, y=[-750]*len(tmp), marker='|', c='black', s=200) #plot lick raster
ax[0].axhline(-750, color='grey', linewidth=1, alpha=0.5, zorder=0)
ax[1].plot(x_times, scnd.gaussian_filter1d(hpc_lfp, sigma=20), color=plt_color, linewidth=4) #plot gaussian filtered V1 vep
ax[1].scatter(x=tmp, y=[-750]*len(tmp), marker='|', c='black', s=200) #plot lick raster
ax[1].axhline(-750, color='grey', linewidth=1, alpha=0.5, zorder=0)
ax[0].set_xlabel('Time (s)')
ax[1].set_xlabel('Time (s)')
ax[0].set_ylabel('V1')
ax[0].set_ylim([-900,(rast_start+unit_gap*5+400)])
ax[0].set_yticks([-500,0,500,rast_start,(rast_start+unit_gap*5)])
ax[0].set_yticklabels([-500,0,500,0,5])
# ax[1].set_ylim([-900,900])
ax[0].set_title(f"{et} ~~ Trial {plt_tr} ~~ {plt_stim} ~~ {tr_beh}")
ax[1].set_ylabel('HPC')
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#plot the unit rasters (~6 units per region - some don't fire on every trial)
unit_spk_tr = unit_spk[unit_spk.trial==plt_tr]
my_pal = sns.color_palette("tab10")
for idx in range(len(v1_units)):
ax[0].axhline(rast_start+unit_gap*idx, color='grey', linewidth=1, alpha=0.5, zorder=0)
for idx in range(len(hpc_units)):
ax[1].axhline(rast_start+unit_gap*idx, color='grey', linewidth=1, alpha=0.5, zorder=0)
for u,uu in unit_spk_tr.groupby('cuid'):
if uu.region.unique()[0]=='v1':
ax[0].scatter(x=uu.trial_spikes, y=[rast_start+unit_gap*v1_units.index(u)]*len(uu.trial_spikes), marker='|', s=120, color=my_pal[v1_units.index(u)])
elif uu.region.unique()[0]=='hippo':
ax[1].scatter(x=uu.trial_spikes, y=[rast_start+unit_gap*hpc_units.index(u)]*len(uu.trial_spikes), marker='|', s=120, color=my_pal[hpc_units.index(u)])
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
sns.despine()
if is_saving:
plt.savefig(rf"C:\Users\AChub_Lab\Desktop\oneEToneTrial_LFPUnitLick_{dd.group.unique()[0]}_tr{plt_tr}_{plt_stim}_{tr_beh}.pdf", transparent=True)
plt.show()
break
# %% [markdown]
# # ! Spike-Phase coherence
# %%
print("LFP dataFrame")
display(all_trial_lfp.head(2))
print("LFP dataFrame to map behavior responses")
display(lfp_fft.head(2))
print("Unit spike time dataFrame")
display(unit_spikes.head(2))
print("Unit psth for visually responsive units")
display(vr_units.head(2))
# %%
def butter_bandpass_filter(mydata, lowcut, highcut, fs, order=4):
nyq = 0.5 * fs
low = float(lowcut)/nyq
high = float(highcut)/nyq
b, a = ssig.butter(order, [low,high], btype='band')
y = ssig.filtfilt(b, a, mydata)
return y
# %%
#This cell cuts the unit spike time DataFrame to only include visually responsive units from the vr_psth DataFrame
vr_unit_ls = vr_units.cuid.unique()
vr_unit_spikes = unit_spikes[unit_spikes['cuid'].isin(vr_unit_ls)]
# %%
#### THIS SPC LOOKS AT LFP CHANNEL THATS SPECIFIC TO THE UNIT DEPTH
lfp_sr = 2500
unit_sr = 30000
time_win = [0.5, 1.5]
# v1_ch = 260
# hpc_ch = 110
spk_ph_coh_ls, et_num = [], 1
for d,dd in all_trial_lfp.groupby(["et","trial"]):
#progress bar for looping through ets
if (d[1]+1)%150==0:
print(f"Done with {et_num} out of {all_trial_lfp.et.nunique()} mice")
et_num+=1
#metadata for saving below
et_group = dd.group.unique()[0]
tr_stim = 'Go+' if dd.stim_id.unique()[0]=='0' else ('Go-' if dd.stim_id.unique()[0]=='1' else 'No-Go')
et_behav = lfp_fft[(lfp_fft.et==d[0])&(lfp_fft.trial==d[1])].behav.unique()[0] #this is the et:trial behavior response
# load spikes df, and limit to the specific et
et_spks = vr_unit_spikes[(vr_unit_spikes.et==d[0])&(vr_unit_spikes.trial==d[1])&(vr_unit_spikes.trial_spikes>=time_win[0])&(vr_unit_spikes.trial_spikes<=time_win[1])] #trim units df to time window
et_spks_v1hpc = et_spks[(et_spks.region=='v1')|(et_spks.region=='hippo')]
for u,uu in et_spks_v1hpc.groupby('cuid'):
unit_depth_ch = int(384-(uu.depth.unique()[0])/10)
my_region = 'V1' if uu.region.unique()[0] == 'v1' else 'HPC'
tr_lfp = dd.lfp_data.values[0][unit_depth_ch] #pick LFP channel that relates to the unit depth
filt_lfp = butter_bandpass_filter(tr_lfp, lowcut=4, highcut=8, fs=lfp_sr, order=4) #filter to theta range
reg_lfp_phase = np.angle(ssig.hilbert(filt_lfp)) #convert to phase domain [-pi, pi]
u_spk = uu.trial_spikes.values
unit_phase_ls = []
my_time = np.linspace(time_win[0],time_win[1],2500)
for i in range(len(my_time)-1):
for spk in u_spk:
if (spk<my_time[i+1])&(spk>=my_time[i]):
unit_phase_ls.append(reg_lfp_phase[i])
u_phs_circm = sstat.circmean(unit_phase_ls, high=np.pi, low=-np.pi)
spk_ph_coh_ls.append(pd.DataFrame({'trial':[d[1]], 'phase':[u_phs_circm], 'cuid':[u], 'et':[d[0]],
'group':[et_group], 'region':[my_region], 'stim_id':[tr_stim], 'behav':[et_behav]}))
spk_ph_coh = pd.concat(spk_ph_coh_ls, ignore_index=True)
spk_ph_coh.head()
# %%
#takes the circular mean of each unit, split by stim_id and behavior
spk_ph_coh_mean = []
for d,dd in spk_ph_coh.groupby(['cuid','stim_id','behav']):
mean_phs_circm = sstat.circmean(dd.phase.values, high=np.pi, low=-np.pi)
spk_ph_coh_mean.append(pd.DataFrame({'mean_phase':[mean_phs_circm], 'cuid':[d[0]], 'et':dd.et.unique()[0],
'group':dd.group.unique()[0], 'region':dd.region.unique()[0], 'stim_id':[d[1]], 'behav':[d[2]]}))
spk_ph_coh_mean = pd.concat(spk_ph_coh_mean, ignore_index=True)
spk_ph_coh_mean.head()
# %%
for d,dd in spk_ph_coh_mean.groupby('stim_id'):
g=sns.displot(data=dd, x='mean_phase', col='behav', col_order=['CR', 'FA'] if d=='No-Go' else ['H', 'M'], row='region', row_order=['V1', 'HPC'],
hue='group', hue_order=['WT','FX'], palette={'WT':'cyan', 'FX':'magenta'},
kde=True, rug=False, stat='density', binwidth=(np.pi/10), height=4, aspect=1.25, linewidth=0)
plt.suptitle(d)
g.set(xticks=[-np.pi, -np.pi/2, 0, np.pi/2, np.pi], yticks=[0, 0.03, 0.06], ylim=[0,0.065])
g.set_xticklabels([r'- $\pi$',r'- $\frac{\pi}{2}$', '0', r'$\frac{\pi}{2}$', r'$\pi$'])
print("2-sided KS test comparing the distributions of WT & FX")
for e,ee in dd.groupby(['region','behav']):
x,y = ee[ee.group=='WT']['mean_phase'].values, ee[ee.group=='FX']['mean_phase'].values # N = #units
res = sstat.ks_2samp(x,y)
plt_stats = "***" if res.pvalue < 0.001 else ("**" if res.pvalue<0.01 else ("*" if res.pvalue <0.05 else "ns"))
print(f"{d} {e} -- N: WT={len(x)}, FX={len(y)} -- p={res.pvalue} --", plt_stats)
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\Spike_phase_coherence\SPC_GenotypeSplit_histogram_{d}.pdf", transparent=True)
plt.show()
# %%
for d,dd in spk_ph_coh_mean.groupby('stim_id'):
g=sns.displot(data=dd, x='mean_phase', col='group', col_order=['WT','FX'], row='region', row_order=['V1', 'HPC'],
hue='behav', palette={'H':'blue', 'M':'grey', 'CR':'blue', 'FA':'grey'},
kde=True, rug=False, stat='density', binwidth=(np.pi/10), height=4, aspect=1.25, linewidth=0)
plt.suptitle(d)
g.set(xticks=[-np.pi, -np.pi/2, 0, np.pi/2, np.pi], yticks=[0, 0.03, 0.06], ylim=[0,0.065])
g.set_xticklabels([r'- $\frac{\pi}{2}$', r'- $\pi$', '0', r'$\frac{\pi}{2}$', r'$\pi$'])
print("2-sided KS test comparing the distributions of behav within genotype")
for e,ee in dd.groupby(['region','group']):
if d == 'No-Go':
x,y = ee[ee.behav=='CR']['mean_phase'].values, ee[ee.behav=='FA']['mean_phase'].values # N = #units
else:
x,y = ee[ee.behav=='H']['mean_phase'].values, ee[ee.behav=='M']['mean_phase'].values # N = #units
res = sstat.ks_2samp(x,y)
plt_stats = "***" if res.pvalue < 0.001 else ("**" if res.pvalue<0.01 else ("*" if res.pvalue <0.05 else "ns"))
print(f"{d} {e} -- N: H/CR={len(x)}, M/FA={len(y)} -- p={res.pvalue} --", plt_stats)
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\temp_figures\Spike_phase_coherence\SPC_BehaviorSplit_histogram_{d}.pdf", transparent=True)
plt.show()
# %% [markdown]
# # ! Aside, more unit analysis
# This is an aside for the above cells. I'm working on making a 2d plot for the unit responses split by behavior
# %% [markdown]
# ## PSTH dataframe creation
# %%
print("LFP dataFrame to map behavior responses")
display(lfp_fft.head(2))
print("Unit spike time dataFrame")
display(unit_spikes.head(2))
print("Unit psth (no behavior) used to id visually responsize units")
display(vr_units.head(2))
# %%
et_trl_behav_dict = {}
for t,tt in lfp_fft.groupby('et'):
inner = {}
for m,mm in tt.groupby('trial'):
inner.update({m:mm.behav.unique()[0]})
et_trl_behav_dict.update({t:inner})
# print(et_trl_behav_dict)
# %%
for t,tt in enumerate(unit_spikes.et.unique()):
unit_spikes.loc[(unit_spikes.et==tt), 'behav'] = unit_spikes[unit_spikes.et==tt].trial.map(et_trl_behav_dict[tt])
unit_spikes.head()
# %%
#This cell cuts the unit spike time DataFrame to only include visually responsive units from the vr_psth DataFrame
#also only includes units that are in V1 or HPC
vr_unit_ls = vr_units.cuid.unique()
region_ls = ['v1','hippo']
vr_unit_spikes = unit_spikes[(unit_spikes['cuid'].isin(vr_unit_ls))&(unit_spikes['region'].isin(region_ls))]
# %%
ls_psth = []
th_bin = 0.01
trial_length = 3.0
num_units, un_idx = vr_unit_spikes['cuid'].nunique(), 0
for l,ll in vr_unit_spikes.groupby(['cuid', 'stim_id', 'behav']): ##### I changed this from df_rez to data_df to check the units
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
trials_number_not_empty = len(ll.trial.unique())
h, ttr = mz_ena.PSTH(ll.trial_spikes, th_bin, trial_length, trials_number_not_empty)
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
zscore = sstat.mstats.zscore(h)
mean = np.mean(h[0:50])#The Baseline period. Be sure it matches time course of experiments##
std = 1 if mean<=0 else np.std(h[0:50])
ztc = (h - mean)/std
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
tr_stim = 'Go+' if l[1]==0 else ('Go-' if l[1]==1 else 'No-Go')
my_region = 'V1' if ll.region.unique()[0] == 'v1' else 'HPC'
tr_group = 'WT' if ll.group.unique()[0] == "A" else "FX"
df_psth_tmp = pd.DataFrame({'times':ttr, 'stim_id':tr_stim, 'Hz':h, 'depth':ll.depth.unique()[0],
'zscore':zscore, 'ztc':ztc, 'et':ll.et.unique()[0], 'cc': ll.cc.unique()[0],
'cuid':l[0], 'region':my_region, 'group':tr_group, 'behav':l[2]})
ls_psth.append(df_psth_tmp)
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
behav_psth = pd.concat(ls_psth)
behav_psth.head()
# %% [markdown]
# ## Power Spectrum Analysis, unit activity, split by behavior
# %%
def get_unit_fft(unit_df):
unit_arr = unit_df[(unit_df.times>0.5)&(unit_df.times<1.8)].zscore.values
freq = np.arange(unit_arr.shape[0]) / unit_arr.shape[0] * 100
freq = freq[:freq.shape[0]//2]
f = np.fft.fft(unit_arr)
magnitude_spectrum = (np.abs(f)[:freq.shape[0]])
return freq, magnitude_spectrum
# %%
behav_units_fft = []
for unit,df in behav_psth.groupby(['stim_id','behav','cuid']):
fft_freq, fft_amp = get_unit_fft(df)
bands = [[4,8]]
auc_val,band_val = [],[]
for ran in bands:
lower = fft_freq.searchsorted(ran[0], 'left')
upper = fft_freq.searchsorted(ran[1], 'right') -1
val = auc(fft_freq[lower:upper],fft_amp[lower:upper])
auc_val.append(val)
band_val.append(str(ran))
behav_units_fft.append(pd.DataFrame({'band': band_val, 'auc': auc_val, 'cuid': unit[-1],
'stim_id': unit[0], 'group': df.group.unique()[0], 'behav':unit[1],
'region': df.region.unique()[0], 'et': df.et.unique()[0]}))
behav_units_fft = pd.concat(behav_units_fft, ignore_index=True)
behav_units_fft.head()
# %%
# this is finding the threshold for each group/stim_id pairing to threshold units in cell below
num_std = 1.0
nested_dict={}
for d,dd in behav_units_fft.groupby('group'):
nested_dict[d] = {}
for stim in ["Go", "No-Go"]:
ee = dd[dd.stim_id!="No-Go"] if stim=="Go" else dd[dd.stim_id=="No-Go"]
thresh = ee.auc.values.mean()+num_std*ee.auc.values.std() # threshold of mean + 1 std.
nested_dict[d][stim] = thresh
nested_dict
# %%
unit_thresh_dict = {}
for stim in ["Go", "No-Go"]:
temp_df = behav_units_fft[behav_units_fft.stim_id!="No-Go"] if stim=="Go" else behav_units_fft[behav_units_fft.stim_id=="No-Go"]
unit_thresh_dict[stim] = {}
for d,dd in temp_df.groupby(['group','cuid']):
thresh = nested_dict[d[0]][stim]
unit_thresh_dict[stim][d[1]] = 'yes' if dd.auc.values.mean() >= thresh else 'no'
# print(unit_thresh_dict)
# %%
# Function to categorize values
def categorize(value):
if value == 'No-Go':
return 'No-Go'
else:
return 'Go'
behav_units_fft['stim_id2'] = behav_units_fft['stim_id'].apply(categorize)
behav_units_fft['thresh'] = behav_units_fft.apply(lambda x: unit_thresh_dict[x['stim_id2']][x['cuid']], axis=1) #map the nexted dictionary to the stim/cuid pairing
behav_units_fft.head()
# %%
for d,dd in behav_units_fft[behav_units_fft.thresh=='yes'].groupby('stim_id2'):
print("~~~~~~~~~~ comparing behavior conditions within each group ~~~~~~~~~~")
for e,ee in dd.groupby(['region', 'group']):
conds = ee.behav.unique()
x=ee[ee.behav==conds[0]].auc.values
y=ee[ee.behav==conds[1]].auc.values
res = sstat.mannwhitneyu(x, y)
pstar = '***' if res.pvalue < 0.001 else ('**' if res.pvalue < 0.01 else ('*' if res.pvalue < 0.05 else "ns"))
print(f"{e} -- {res} -- {pstar}")
print("~~~~~~~~~~ comparing groups within each behavior condition ~~~~~~~~~~")
for e,ee in dd.groupby(['region', 'behav']):
conds = ee.group.unique()
x=ee[ee.group==conds[0]].auc.values
y=ee[ee.group==conds[1]].auc.values
res = sstat.mannwhitneyu(x, y)
pstar = '***' if res.pvalue < 0.001 else ('**' if res.pvalue < 0.01 else ('*' if res.pvalue < 0.05 else "ns"))
print(f"{e} -- {res} -- {pstar}")
g = sns.catplot(dd, kind="bar",
x="behav", y="auc", hue='group', hue_order=['WT','FX'],
col='region', col_order=['V1', 'HPC'],
height=4, aspect=1.2)
g.set_ylabels(f"{d} AUC")
plt.show()
# %%
# %%
# %% [markdown]
# ## Plot it
# %%
osc_characteristic = 'auc'
unit_evidence = []
for d,dd in behav_units_fft[(behav_units_fft.et.isin(active_ets))].groupby(['et','region']):
go_c = dd[(dd.stim_id=='Go+')&(dd.behav=='H')|(dd.stim_id=='Go-')&(dd.behav=='H')][osc_characteristic].mean()
go_i = dd[(dd.stim_id=='Go+')&(dd.behav=='M')|(dd.stim_id=='Go-')&(dd.behav=='M')][osc_characteristic].mean()
go_all = dd[(dd.stim_id=='Go+')|(dd.stim_id=='Go-')][osc_characteristic].mean()
ng_c = dd[(dd.stim_id=='No-Go')&(dd.behav=='CR')][osc_characteristic].mean()
ng_i = dd[(dd.stim_id=='No-Go')&(dd.behav=='FA')][osc_characteristic].mean()
ng_all = dd[dd.stim_id=='No-Go'][osc_characteristic].mean()
go = (go_c-go_i)/(go_c+go_i)
ng = (ng_c-ng_i)/(ng_c+ng_i)
g_ng = (go-ng)/(go+ng)
go_base = (go_c-go_all)/(dd[(dd.stim_id=='Go+')|(dd.stim_id=='Go-')][osc_characteristic].std())
ng_base = (ng_c-ng_all)/(dd[dd.stim_id=='No-Go'][osc_characteristic].std())
foo = lfp_fft[(lfp_fft.et==d[0])&(lfp_fft.region==d[-1])]
lick_dict = foo.groupby('behav').trial.nunique().to_dict()
d_prime = sstat.norm.ppf(lick_dict['H']/100) - sstat.norm.ppf(lick_dict['FA']/50)
unit_evidence.append(pd.DataFrame({'goC':go_c, 'goIC':go_i,'ngC':ng_c, 'ngIC':ng_i,'go_avg':go, 'ng_avg':ng, 'gng_avg':g_ng, 'go_base':go_base, 'ng_base':ng_base,
'go_licks':lick_dict['H']/100, 'ng_licks':lick_dict['CR']/50, 'gng_licks':(lick_dict['H']+lick_dict['CR'])/150, 'Dprime':d_prime,
'go2_licks':lick_dict['M']/100, 'ng2_licks':lick_dict['FA']/50,
'et':d[0], 'group':dd.group.unique()[0], 'region':d[-1]}, index=[0]))
unit_evidence = pd.concat(unit_evidence, ignore_index=True)
unit_evidence.head()
# %%
dot_size = 50 # 50, aesthetic choice for the size of the scatter plot dots
plot_linReg = False # adds a linear regression to the 2d scatter
plot_ell = True #adds the covariance confidence ellipse to the data
num_std = 1.5 # standard deviation you want to have contained within the ellipse
# for val in [['go_licks', 'go_base'], ['ng_licks', 'ng_base']]:
# for val in [['go_licks', 'go_avg'], ['ng_licks', 'ng_avg']]:#, ['gng_licks', 'gng_avg'], ['Dprime', 'gng_avg']]:
for val in [['go_licks', 'goC'], ['ng_licks', 'ngC']]:
# for val in [['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]:
my_x, my_y = val[0], val[1]
fig, ax = plt.subplots(1, 2, figsize=(10, 4), sharex=True, sharey=True)
for idx, reg in enumerate(['V1','HPC']):
x = unit_evidence[(unit_evidence.group=='WT')&(unit_evidence.region==reg)][my_x].values
y = unit_evidence[(unit_evidence.group=='WT')&(unit_evidence.region==reg)][my_y].values
x1 = unit_evidence[(unit_evidence.group=='FX')&(unit_evidence.region==reg)][my_x].values
y1 = unit_evidence[(unit_evidence.group=='FX')&(unit_evidence.region==reg)][my_y].values
# plot linear regression of the points & print R**2 value of the fit
if plot_linReg:
res, res1 = sstat.linregress(x, y), sstat.linregress(x1, y1)
print(f"------------- {my_x} -------------")
print(f'{reg} - Slope: WT {res.slope:.6f}, FX {res1.slope:.6f}')
print(f'{reg} - R squared: WT {res.rvalue**2:.6f}, FX {res1.rvalue**2:.6f}')
#"The p-value for a hypothesis test whose null hypothesis is that the slope is zero, using Wald Test with t-distribution of the test statistic."
print(f'{reg} - pval: WT {res.pvalue:.6f}, FX {res1.pvalue:.6f}')
ax[idx].plot(x, res.intercept+res.slope*x, color='cyan')
ax[idx].plot(x1, res1.intercept+res1.slope*x1, color='magenta')
# Plot scatter
ax[idx].scatter(x, y, s=dot_size, edgecolors=None, color='cyan', label='WT')
ax[idx].scatter(x1, y1, s=dot_size, edgecolors=None, color='magenta', label='FX')
# Plot ellipse using above function
if plot_ell:
confidence_ellipse(x, y, ax[idx], n_std=num_std, edgecolor=None, facecolor='cyan', alpha=0.3)
confidence_ellipse(x1, y1, ax[idx], n_std=num_std, edgecolor=None, facecolor='magenta', alpha=0.3)
if my_x == 'gng_licks':
ax[idx].set_xlim([0,1.1])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("Correct %")
ax[0].set_ylabel("GNG Theta")
elif my_x == 'Dprime':
ax[idx].set_xlim([-1,2])
# ax[idx].set_ylim([0.45,0.57])
ax[idx].set_xlabel("d'")
ax[0].set_ylabel("GNG Theta")
else:
ax[idx].set_xlim([-0.1,1.1])
ax[idx].set_xlabel("Correct rate")
if (my_y=='go_avg')|(my_y=='ng_avg'):
ax[0].set_ylabel("Go Theta") if my_x=='go_licks' else ax[0].set_ylabel("No-Go Theta")
else:
ax[0].set_ylabel("Theta power")
if (my_y=='goIC')|(my_y=='ngIC'):
ax[idx].set_ylim([16,39])
ax[idx].set_xlabel("Incorrect rate")
ax[idx].set_title(f"{reg} units")
plt.tight_layout()
sns.despine()
# plt.savefig(rf"C:\Users\AChub_Lab\Desktop\2dscatter_units_IncorrectTrials_{my_x}_v1hpc.pdf", transparent=True)
plt.show()
# %%
# Calculate a Spearman correlation coefficient with associated p-value on the 2d scatter plots
def statistic(x): # permute only `x`
return sstat.spearmanr(x, y).statistic
# for val in [['go_licks', 'go_avg'], ['ng_licks', 'ng_avg'], ['gng_licks', 'gng_avg'], ['Dprime', 'gng_avg']]:
for val in [['go_licks', 'goC'], ['ng_licks', 'ngC'], ['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]: #correct/incorrect trial responses only
my_x, my_y = val[0], val[1]
for idx, reg in enumerate(['V1','HPC']):
for g in ['WT', 'FX']:
x = unit_evidence[(unit_evidence.group==g)&(unit_evidence.region==reg)][my_x].values
y = unit_evidence[(unit_evidence.group==g)&(unit_evidence.region==reg)][my_y].values
res_exact = sstat.permutation_test((x,), statistic, permutation_type='pairings')
print(val, reg, g)
pval = "***" if res_exact.pvalue < 0.001 else ("**" if res_exact.pvalue<0.01 else ("*" if res_exact.pvalue <0.05 else "ns"))
print(f"stat = {res_exact.statistic:.6f} -- pval = {res_exact.pvalue:.6f} -- {pval}")
# %%
# trying a 2-dimensional, 2-sample KS test for comparing the WT and FX groups (https://github.com/syrte/ndtest)
import ndtest
for val in [['go_licks', 'goC'], ['ng_licks', 'ngC'], ['go2_licks', 'goIC'], ['ng2_licks', 'ngIC']]: #correct/incorrect trial responses only
my_x, my_y = val[0], val[1]
for idx, reg in enumerate(['V1','HPC']):
x = unit_evidence[(unit_evidence.group=='WT')&(unit_evidence.region==reg)][my_x].values
y = unit_evidence[(unit_evidence.group=='WT')&(unit_evidence.region==reg)][my_y].values
x1 = unit_evidence[(unit_evidence.group=='FX')&(unit_evidence.region==reg)][my_x].values
y1 = unit_evidence[(unit_evidence.group=='FX')&(unit_evidence.region==reg)][my_y].values
P, D = ndtest.ks2d2s(x, y, x1, y1, extra=True)
pval = "***" if P<0.001 else ("**" if P<0.01 else ("*" if P<0.05 else "ns"))
print(f"{my_y} -- {reg} -- d={D:.7f} -- p={P:.7f} -- {pval}")
# %%

Behavior_Theta_Correlation.ipynb at commit 05fe09d, no license · at the source

Overview

Authors: Michael P Zimmerman^1,2, Mowen Yin^1,2, Kevin R Cragg^1,2, Sanghamitra Nareddula¹, Varun M Kumar³, Violeta Saldarriaga¹, Adriana Rotger², Rachel Lehman¹, Paige Edens¹, Caroline Powell¹, Jenna Barry¹, Sein Kim¹, Joseph G Makin³, Alexander A Chubykin^1,2,4

ORCID iDs: Mowen Yin, Paige Edens, Alexander A Chubykin

Institutions: Purdue University West Lafayette (United States)

Journal: Cell reports, volume 45, issue 7, article 117590

Dates: published online 20 June 2026; in print 28 July 2026

Type: Research article · Language: English

License: CC BY-NC

Identifiers: DOI 10.1016/j.celrep.2026.117590 · PMID 42322609 · PMCID PMC13472280 · OpenAlex W7165403289

Open access: gold, a free copy (OpenAlex)

Status: code verified

Categories: EEG (modality), mouse (organism), other condition (population), autism (population), systems (subfield)

Methods: Spectral & time-frequency, Preprocessing, Statistics, Machine learning, Connectivity, Evoked potentials, Single-unit activity, calcium imaging

Keywords: Memory, Learning, Hippocampus, Vision, Primary visual cortex, prefrontal cortex, Autism, Fragile X Syndrome, Neural Circuits, Cp: Neuroscience

MeSH: Behavior, Animal*, Fragile X Messenger Ribonucleoprotein 1*, Fragile X Syndrome*, Theta Rhythm*, Visual Cortex*, Animals, Hippocampus, Male, Mice, Mice, Inbred C57BL, Mice, Knockout, Prefrontal Cortex (* major topic)

Topic: Memory and Neural Mechanisms (Cognitive Neuroscience, Neuroscience), according to OpenAlex

Funding: NIMH (R01 MH116500); National Institutes of Health; NIMH NIH HHS (R01 MH116500)

Citations: not cited yet (Europe PMC); 105 references in the paper

Research resources: Alexa Anti-Chicken 488 RRID:AB_142924, Wild type C57BL/6 mice RRID:IMSR_JAX:000664, B6.129P2-Fmr1tm1Cgr/J (Fmr1 KO mice) RRID:IMSR_JAX:003024, RRID:IMSR_JAX:008069, RRID:IMSR_JAX:012569

Abstract

Fragile X syndrome (FX) is associated with sensory processing and learning deficits. Visual familiarity evokes persistent theta oscillations during passive stimulus presentation in the primary visual cortex (V1) and the hippocampus (HPC), which are impaired in V1 of FX. How does this activity change during active behavior? To address this, we performed Neuropixels recordings in V1, HPC, and the prefrontal cortex (PFC) during Go/No-Go visual discrimination behavior. Theta oscillations are reduced in FX in both V1 and HPC and are abolished during No-Go trials, correlating with excessive, incorrect licking behavior. In wild-type (WT) mice, V1 theta power strongly correlates with correct behavioral outcomes. PFC shows significantly reduced cue-related responses in FX. Together, these findings show loss of behavioral inhibition in FX correlated with attenuated theta activity in V1 and HPC and deficient top-down control from PFC. This work sheds light on circuit-level impairments underlying behavioral deficits in FX for potential therapeutic interventions.

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

Repositories

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

chubykin/V1-HPC-PFC-Behavior-Oscillations

License: none: the authors keep all their rights

State: the link answers, verified on 27 September 2026

Evidence: files inventoried

Commit: 05fe09d84329d3582aee3ca74492e9b623a33d74, 19 May 2026

Languages: Jupyter (44), Python (24)

Size: 69 files, 68 scripts

Software Heritage: not archived

Found in: “Data and code availability”

Holds: README, 44 notebooks

Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation

Tools: NumPy (29 files), SciPy (27 files), pandas (21 files), Matplotlib (20 files), seaborn (15 files), statsmodels (12 files), Pingouin (4 files), scikit-learn (4 files), Neo (2 files), Open Ephys analysis tools (1 file)

Availability: 1 check, the latest on 27 September 2026: the link answers

27 September 2026: the link answers

29 files

Behavior_Theta_Correlation/Behavior_Theta_Correlation.ipynb, Jupyter, 1,672 lines, 8 matches
LFP/64Chs_LFP/64Ch_lfp.ipynb, Jupyter, 1,107 lines, 2 matches
LFP/Neuropixels_LFP/1_saving_LFP_arrays.ipynb, Jupyter, 307 lines
LFP/Neuropixels_LFP/2_combo_intervals.ipynb, Jupyter, 273 lines
LFP/Neuropixels_LFP/2_reward_LFP-novel.ipynb, Jupyter, 430 lines, 2 matches
LFP/Neuropixels_LFP/2_reward_LFP.ipynb, Jupyter, 695 lines
LFP/Neuropixels_LFP/3_TF_analysis-novel.ipynb, Jupyter, 293 lines, 1 match
LFP/Neuropixels_LFP/3_TF_analysis.ipynb, Jupyter, 435 lines
LFP/Neuropixels_LFP/4_CSD_analysis.ipynb, Jupyter, 664 lines
LFP/Neuropixels_LFP/4_CSD_analysis_novel.ipynb, Jupyter, 403 lines
LFP/Neuropixels_LFP/Python3_OpenEphys_V14.py, Python, 461 lines
LFP/Neuropixels_LFP/Python3_OpenEphys_orig.py, Python, 731 lines
LFP/Neuropixels_LFP/Python3_OpenOE_AC_map_functions_v1_08_30s.py, Python, 1,672 lines
LFP/Neuropixels_LFP/Python3_icsd.py, Python, 930 lines
LFP/Neuropixels_LFP/Python3_new_OpenEphys_orig.py, Python, 461 lines
LFP/Neuropixels_LFP/mz_LFP_functions.py, Python, 390 lines
LFP/V1_HPC/1_HPC_LFP.ipynb, Jupyter, 397 lines, 1 match
LFP/V1_HPC/2_HPC_TF.ipynb, Jupyter, 392 lines, 1 match
LFP/V1_HPC/3_TrlByTrl_analysis.ipynb, Jupyter, 638 lines
LFP/V1_HPC/New_HPC_LFP_traces.ipynb, Jupyter, 187 lines
LFP/V1_HPC/Python3_OpenEphys_V14.py, Python, 461 lines
LFP/V1_HPC/Python3_OpenEphys_orig.py, Python, 731 lines
LFP/V1_HPC/Python3_OpenOE_AC_map_functions_v1_08_30s.py, Python, 1,672 lines
LFP/V1_HPC/Python3_icsd.py, Python, 930 lines
LFP/V1_HPC/Python3_new_OpenEphys_orig.py, Python, 461 lines
LFP/V1_HPC/mz_LFP_functions.py, Python, 390 lines, 1 match
LFP/V1_HPC/plv_phd_afterSaving.ipynb, Jupyter, 683 lines, 2 matches
ML Code/Create_spikecount_df.ipynb, Jupyter, 81 lines
ML Code/accvsunits.py, Python, 184 lines, 2 matches
repository limit reached (2,000 files or 30 MB): the rest is at the source (40 files)

At the source: github.com/chubykin/V1-HPC-PFC-Behavior-Oscillations

cortex-lab/allenCCF

License: none: the authors keep all their rights

State: the link answers, verified on 27 September 2026

Evidence: files inventoried

Commit: e5a57fe7e1c9fb333fec51c29a8471131c233a76, 15 July 2025

Languages: MATLAB (61)

Size: 75 files, 61 scripts

Software Heritage: not archived

Found in: the resources table

Holds: README

Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation

Tools: Image Processing Toolbox (7 files), Statistics and Machine Learning Toolbox (1 file)

Availability: 1 check, the latest on 27 September 2026: the link answers

27 September 2026: the link answers

62 files

Browsing Functions/AtlasTransformBrowser.m, MATLAB, 1,196 lines
Browsing Functions/CCF_to_FP.m, MATLAB, 66 lines
Browsing Functions/addAllenCtxOutlines.m, MATLAB, 51 lines
Browsing Functions/aggregateAcr.m, MATLAB, 146 lines
Browsing Functions/allenAtlasBrowser.m, MATLAB, 658 lines
Browsing Functions/allenAtlasBrowser_original.m, MATLAB, 262 lines
Browsing Functions/allenCCFbregma.m, MATLAB, 8 lines
Browsing Functions/allenTilt.m, MATLAB, 172 lines
Browsing Functions/allen_ccf_2pi.m, MATLAB, 780 lines
Browsing Functions/allen_ccf_colormap.m, MATLAB, 10 lines
Browsing Functions/allen_ccf_npx.m, MATLAB, 1,087 lines
Browsing Functions/allen_ccf_npx_4shank.m, MATLAB, 931 lines
Browsing Functions/allen_ccf_npx_4shank_spherical.m, MATLAB, 987 lines
Browsing Functions/best_fit_line.m, MATLAB, 28 lines
Browsing Functions/customVectorSlice.m, MATLAB, 38 lines
Browsing Functions/distinguishable_colors.m, MATLAB, 152 lines
Browsing Functions/get_offset_map.m, MATLAB, 41 lines
Browsing Functions/gridIn3D.m, MATLAB, 80 lines
Browsing Functions/hierarchicalSelect.m, MATLAB, 130 lines
Browsing Functions/idRegionByAcr.m, MATLAB, 27 lines
Browsing Functions/isAreaOrContains.m, MATLAB, 31 lines
Browsing Functions/loadCCFtoFP.m, MATLAB, 17 lines
Browsing Functions/loadFPtable.m, MATLAB, 16 lines
Browsing Functions/loadStructureTree.m, MATLAB, 48 lines
Browsing Functions/makeSTtree.m, MATLAB, 45 lines
Browsing Functions/makeSmoothCoords.m, MATLAB, 23 lines
Browsing Functions/natsort.m, MATLAB, 330 lines
Browsing Functions/natsortfiles.m, MATLAB, 169 lines
Browsing Functions/plotAVoverlay.m, MATLAB, 28 lines
Browsing Functions/plotAVslice.m, MATLAB, 13 lines
Browsing Functions/plotAsProbe.m, MATLAB, 37 lines
Browsing Functions/plotBrainGrid.m, MATLAB, 35 lines
Browsing Functions/plotBrainOutlinesByAxis.m, MATLAB, 23 lines
Browsing Functions/plotDistToNearest.m, MATLAB, 123 lines
Browsing Functions/plotDistToNearestToTip.m, MATLAB, 252 lines
Browsing Functions/plotLabelsAsProbe.m, MATLAB, 110 lines
Browsing Functions/plotNeuronOnSliceFromCoord.m, MATLAB, 20 lines
Browsing Functions/plotTVslice.m, MATLAB, 12 lines
Browsing Functions/plotTopDownOutlines.m, MATLAB, 100 lines
Browsing Functions/sagittalSlices.m, MATLAB, 50 lines
Browsing Functions/sanitizeStructureTree.m, MATLAB, 31 lines
Browsing Functions/script_sliceMovie.m, MATLAB, 66 lines
Browsing Functions/selectStructure.m, MATLAB, 264 lines
Browsing Functions/sliceBrowser.m, MATLAB, 142 lines
Browsing Functions/sliceByVector.m, MATLAB, 55 lines
Browsing Functions/sliceOutlineWithRegion.m, MATLAB, 38 lines
Browsing Functions/sliceOutlineWithRegionVec.m, MATLAB, 74 lines
Browsing Functions/transformed_sliceBrowser.m, MATLAB, 171 lines
Histology Functions/HistologyBrowser.m, MATLAB, 161 lines
Histology Functions/HistologyCropper.m, MATLAB, 105 lines
Histology Functions/SliceFlipper.m, MATLAB, 177 lines
Histology Functions/natsort.m, MATLAB, 330 lines
Histology Functions/natsortfiles.m, MATLAB, 169 lines
SHARP-Track/Analyze_Clicked_Points.m, MATLAB, 131 lines
SHARP-Track/Analyze_ROIs.m, MATLAB, 158 lines
SHARP-Track/Convert_CCF_Coords_to_FP_Regions.m, MATLAB, 117 lines
SHARP-Track/Convert_Clicked_Points_to_FP_coords.m, MATLAB, 145 lines
SHARP-Track/Display_Probe_Track.m, MATLAB, 229 lines
SHARP-Track/Navigate_Atlas_and_Register_Slices.m, MATLAB, 66 lines
SHARP-Track/Process_Histology.m, MATLAB, 133 lines
setup_utils.m, MATLAB, 47 lines
README.md, Text, 72 lines

At the source: github.com/cortex-lab/allenCCF

The paper's code and data availability statement is in the Data section.

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;
90 scripts, each with its path and the digest of its content;
20 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

zenodo:20273683, at Zenodo; found in “Data and code availability”

Data and code availability

All original datasets have been deposited at Zenodo at https://doi.org/10.5281/zenodo.20273683 and are publicly available as of the date of publication.

All original code has been deposited at GitHub and is publicly available at https://github.com/chubykin/V1-HPC-PFC-Behavior-Oscillationshttps://github.com/chubykin/V1-HPC-PFC-Behavior-Oscillations as of the date of publication.

Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.

Reproduced under the paper's license (CC BY-NC), 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 2, 28 September 2026

Publisher: n/a → Cell Press
Authors: added Alexander A Chubykin (0000-0001-8224-9296); removed Alexander A Chubykin

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 14 authors, 10 keywords, 12 MeSH terms, 3 funders, 103 references, 5 RRIDs.

Cite

This paper

Zimmerman, M. P., Yin, M., Cragg, K. R., Nareddula, S., Kumar, V. M., Saldarriaga, V., Rotger, A., Lehman, R., Edens, P., Powell, C., Barry, J., Kim, S., Makin, J. G., & Chubykin, A. A. (2026). Impaired behavioral inhibition in Fmr1 KO mice is linked to disrupted visual cortex theta oscillations. Cell reports, 45(7), 117590. https://doi.org/10.1016/j.celrep.2026.117590

BibTeX

@article{zimmerman2026impaired, author = {Zimmerman, Michael P and Yin, Mowen and Cragg, Kevin R and Nareddula, Sanghamitra and Kumar, Varun M and Saldarriaga, Violeta and Rotger, Adriana and Lehman, Rachel and Edens, Paige and Powell, Caroline and Barry, Jenna and Kim, Sein and Makin, Joseph G and Chubykin, Alexander A}, title = {{Impaired behavioral inhibition in Fmr1 KO mice is linked to disrupted visual cortex theta oscillations}}, journal = {Cell reports}, year = {2026}, month = jun, volume = {45}, number = {7}, pages = {117590}, publisher = {Cell Press}, issn = {2211-1247}, doi = {10.1016/j.celrep.2026.117590}, url = {https://doi.org/10.1016/j.celrep.2026.117590}, pmid = {42322609}, pmcid = {PMC13472280} }

RIS

TY - JOUR AU - Zimmerman, Michael P AU - Yin, Mowen AU - Cragg, Kevin R AU - Nareddula, Sanghamitra AU - Kumar, Varun M AU - Saldarriaga, Violeta AU - Rotger, Adriana AU - Lehman, Rachel AU - Edens, Paige AU - Powell, Caroline AU - Barry, Jenna AU - Kim, Sein AU - Makin, Joseph G AU - Chubykin, Alexander A TI - Impaired behavioral inhibition in Fmr1 KO mice is linked to disrupted visual cortex theta oscillations T2 - Cell reports J2 - Cell Rep PY - 2026 DA - 2026/06/20 VL - 45 IS - 7 SP - 117590 SN - 2211-1247 PB - Cell Press DO - 10.1016/j.celrep.2026.117590 UR - https://doi.org/10.1016/j.celrep.2026.117590 LA - en ER -

CSL-JSON

{ "id": "10.1016/j.celrep.2026.117590", "type": "article-journal", "title": "Impaired behavioral inhibition in Fmr1 KO mice is linked to disrupted visual cortex theta oscillations", "container-title": "Cell reports", "author": [ { "family": "Zimmerman", "given": "Michael P" }, { "family": "Yin", "given": "Mowen" }, { "family": "Cragg", "given": "Kevin R" }, { "family": "Nareddula", "given": "Sanghamitra" }, { "family": "Kumar", "given": "Varun M" }, { "family": "Saldarriaga", "given": "Violeta" }, { "family": "Rotger", "given": "Adriana" }, { "family": "Lehman", "given": "Rachel" }, { "family": "Edens", "given": "Paige" }, { "family": "Powell", "given": "Caroline" }, { "family": "Barry", "given": "Jenna" }, { "family": "Kim", "given": "Sein" }, { "family": "Makin", "given": "Joseph G" }, { "family": "Chubykin", "given": "Alexander A" } ], "container-title-short": "Cell Rep", "volume": "45", "issue": "7", "page": "117590", "DOI": "10.1016/j.celrep.2026.117590", "PMID": "42322609", "PMCID": "PMC13472280", "ISSN": "2211-1247", "publisher": "Cell Press", "URL": "https://doi.org/10.1016/j.celrep.2026.117590", "language": "en", "issued": { "date-parts": [ [ 2026, 6, 20 ] ] } }

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.1038/s41586-026-10501-y [code]: Long-term editing of brain circuits using an engineered electrical synapse.
Journal: Nature
In common: Neo, Pingouin, statsmodels, 7 other tools, mouse, 1 reference
[2] doi:10.1016/j.celrep.2026.117420 [code]: Neural population dynamics of direct electrical stimulation of neocortex.
Journal: Cell reports
In common: Open Ephys analysis tools, statsmodels, seaborn, 5 other tools, systems, mouse, 2 references
[3] doi:10.1016/j.patter.2026.101590 [code]: Density-based longitudinal neuron tracking in high-density electrophysiological recordings.
Journal: Patterns (New York, N.Y.)
In common: Neo, statsmodels, Statistics and Machine Learning Toolbox, 6 other tools, 2 references
[4] doi:10.1038/s41593-026-02314-z [code]: Low-dimensional population dynamics in the brainstem gate REM sleep.
Journal: Nature neuroscience
In common: Pingouin, statsmodels, seaborn, 5 other tools, systems, EEG, mouse, 2 references
[5] doi:10.7554/elife.108408 [code]: Frequency and laminar profile of feature-specific visual activity revealed by interleaved EEG-fMRI.
Journal: eLife
In common: Pingouin, Image Processing Toolbox, statsmodels, 6 other tools, systems, EEG, 1 reference
[6] doi:10.1038/s41593-026-02255-7 [code]: Neural circuits encode prior knowledge of temporal statistics.
Journal: Nature neuroscience
In common: Open Ephys analysis tools, Image Processing Toolbox, Statistics and Machine Learning Toolbox, 5 other tools, systems, mouse, 1 reference
[7] doi:10.1002/aur.70312 [code]: Aberrant Neural Entrainment to Word-Level Speech Patterns in Fragile X Syndrome: Evidence for a Statistical Learning Deficit.
Journal: Autism research : official journal of the International Society for Autism Research
In common: Statistics and Machine Learning Toolbox, autism, EEG, other condition, 5 references
[8] doi:10.1038/s41593-026-02232-0 [code]: Entorhinal cortex represents task-relevant remote locations independently of CA1.
Journal: Nature neuroscience
In common: Pingouin, Image Processing Toolbox, statsmodels, 7 other tools, systems, mouse
[9] doi:10.1002/advs.202519479 [code]: Diminished Signal-to-Noise Ratio Disrupts Somatosensory Population Encoding and Drives Tactile Hyposensitivity in the Fmr1<sup>-/y</sup> Autism Model.
Journal: Advanced science (Weinheim, Baden-Wurttemberg, Germany)
In common: Pingouin, statsmodels, seaborn, 5 other tools, autism, other condition, mouse, 1 reference
[10] doi:10.1126/sciadv.aef0343 [code]: Learning induces activation-mechanism-dependent neural plasticity in an intracortical microstimulation task.
Journal: Science advances
In common: Neo, Image Processing Toolbox, statsmodels, 7 other tools

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.