Skip to content

Commit

Permalink
Initial custom changes commit
Browse files Browse the repository at this point in the history
1. Fernandez-Ruiz 2017 PV vs SST metrics
2. Adjust processing code so doesn't require waveforms to be run
3. Calculate CCG and transmission probability in 20s bins
4. Runners for my project & Kei's paper
  • Loading branch information
emilyasterjones committed Mar 24, 2023
1 parent 271d74a commit 8da690e
Show file tree
Hide file tree
Showing 6 changed files with 166 additions and 12 deletions.
49 changes: 49 additions & 0 deletions +customCalculations/FR2017.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,49 @@
function cell_metrics = FR2017(cell_metrics,session,spikes,parameters)
% This is an example template for creating your own calculations
%
% INPUTS
% cell_metrics - cell_metrics struct
% session - session struct with session-level metadata
% spikes_intervalsExcluded - spikes struct filtered by (manipulation) intervals
% spikes - spikes cell struct
% spikes{1} : all spikes
% spikes{2} : spikes excluding manipulation intervals
% parameters - input parameters to ProcessCellExplorer
%
% OUTPUT
% cell_metrics - updated cell_metrics struct

sr = session.extracellular.sr;
spikes = spikes{1};
if any(spikes.total==0)
cell_indexes = find(spikes.total>0);
else
cell_indexes = 1:spikes.numcells;
end

% acg_mean = nan(1,spikes.numcells);
burst_index = nan(1,spikes.numcells);
refrac = nan(1,spikes.numcells);
for i = cell_indexes
acg = CCG(spikes.times{i},ones(size(spikes.times{i})),'binSize',0.001,'duration',.1,'norm','rate','Fs',1/sr);
acg = acg(51:end); % subset to after spike
% acg_mean(i) = mean(acg);
if max(acg(1:10))==0
burst_index(i) = 0;
else
burst_index(i) = (max(acg(1:10)) - mean(acg(end-10:end)))/max(acg(1:10));
end

[~, max_ind] = max(acg);
for m = 2:max_ind
slope(m) = (acg(m) - acg(1))/(m - 1);
end
slope_sd = std(slope);
exceed_1SD = find(slope > slope_sd);
refrac(i) = exceed_1SD(1) - 1;
end

% cell_metrics.acg_mean_FR = acg_mean;
cell_metrics.burst_index_FR = burst_index;
cell_metrics.refrac_period_FR = refrac;
end
47 changes: 40 additions & 7 deletions ProcessCellMetrics.m
Original file line number Diff line number Diff line change
Expand Up @@ -293,6 +293,11 @@
spikes{1}.cluID = spikes{1}.UID;
end

%%%EAJ added 3/21/2023 in case waveforms not loaded
if ~isfield(spikes{1}, 'shankID')
spikes{1}.shankID = ones(1, length(spikes{1}.UID));
end

if ~isempty(parameters.restrictToIntervals)
if size(parameters.restrictToIntervals,2) ~= 2
error('restrictToIntervals has to be a Nx2 matrix')
Expand Down Expand Up @@ -845,6 +850,28 @@
[trans,prob,prob_uncor,pred] = ce_GetTransProb(ccg2(:,cell_metrics.putativeConnections.inhibitory(i,1),cell_metrics.putativeConnections.inhibitory(i,2)), spikes{spkExclu}.total(cell_metrics.putativeConnections.inhibitory(i,1)), mono_res.binSize, 0.020);
cell_metrics.putativeConnections.inhibitoryTransProb(i) = trans;
end

% EAJ added 3/21/2023: transmission probability binned
cell_metrics.putativeConnections.bins = mono_res.ccgR_bins;
for b = 1:length(mono_res.ccgR_binned)
ccg2 = mono_res.ccgR_binned{b};
ccg2(isnan(ccg2)) = 0;

for i = 1:size(cell_metrics.putativeConnections.excitatory,1)
exc_idx = cell_metrics.putativeConnections.excitatory(i,1);
st_bin_idx = spikes{spkExclu}.times{exc_idx}>mono_res.ccgR_bins(b) & spikes{spkExclu}.times{exc_idx}<=mono_res.ccgR_bins(b+1);
[trans,~,~,~] = ce_GetTransProb(ccg2(:,exc_idx,cell_metrics.putativeConnections.excitatory(i,2)),...
length(spikes{spkExclu}.times{exc_idx}(st_bin_idx)), mono_res.binSize, 0.020);
cell_metrics.putativeConnections.excitatoryTransProb_binned(b,i) = trans;
end
for i = 1:size(cell_metrics.putativeConnections.inhibitory,1)
inh_idx = cell_metrics.putativeConnections.inhibitory(i,1);
st_bin_idx = spikes{spkExclu}.times{inh_idx}>mono_res.ccgR_bins(b) & spikes{spkExclu}.times{inh_idx}<=mono_res.ccgR_bins(b+1);
[trans,~,~,~] = ce_GetTransProb(ccg2(:,inh_idx,cell_metrics.putativeConnections.inhibitory(i,2)),...
length(spikes{spkExclu}.times{inh_idx}(st_bin_idx)), mono_res.binSize, 0.020);
cell_metrics.putativeConnections.inhibitoryTransProb_binned(b,i) = trans;
end
end
else
cell_metrics.putativeConnections.excitatory = [];
cell_metrics.putativeConnections.inhibitory = [];
Expand Down Expand Up @@ -1360,8 +1387,12 @@
cell_metrics.UID = spikes{spkExclu}.UID;
cell_metrics.cluID = spikes{spkExclu}.cluID;
cell_metrics.electrodeGroup = spikes{spkExclu}.shankID;
cell_metrics.maxWaveformCh = spikes{spkExclu}.maxWaveformCh1-1;
cell_metrics.maxWaveformCh1 = spikes{spkExclu}.maxWaveformCh1;
if isfield(spikes{spkExclu},'maxWaveformCh1')
cell_metrics.maxWaveformCh = spikes{spkExclu}.maxWaveformCh1-1;
cell_metrics.maxWaveformCh1 = spikes{spkExclu}.maxWaveformCh1;
else
disp('Waveforms not analyzed; will not be able to determine peak channel')
end
cell_metrics.spikeCount = spikes{spkExclu}.total;

cell_metrics.general.electrodeGroups = session.extracellular.electrodeGroups.channels;
Expand All @@ -1384,7 +1415,7 @@
try
cell_metrics.maxWaveformChannelOrder(j) = find([session.extracellular.electrodeGroups.channels{:}] == spikes{spkExclu}.maxWaveformCh1(j));
catch
warning('peak channel not determined')
%disp('peak channel not determined')
cell_metrics.maxWaveformChannelOrder(j) = nan;
end

Expand Down Expand Up @@ -1438,7 +1469,7 @@
% cell_classification_putativeCellType
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %

if ~isfield(cell_metrics,'putativeCellType') || ~ parameters.keepCellClassification
if isfield(cell_metrics,'putativeCellType') || ~ parameters.keepCellClassification
dispLog(['Cell-type classification schema: ' preferences.putativeCellType.classification_schema],basename)
putativeCellType = celltype_classification.(preferences.putativeCellType.classification_schema)(cell_metrics,preferences);
if all(size(putativeCellType)==[1,cell_metrics.general.cellCount]) && iscell(putativeCellType)
Expand Down Expand Up @@ -1493,9 +1524,11 @@
cell_metrics.general = rmfield(cell_metrics.general,field2remove(test));

% cleaning waveforms
field2remove = {'filt_absolute','filt_zscored', 'raw_absolute','raw_zscored','time_zscored'};
test2 = isfield(cell_metrics.waveforms,field2remove);
cell_metrics.waveforms = rmfield(cell_metrics.waveforms,field2remove(test2));
if isfield(cell_metrics, 'waveforms')
field2remove = {'filt_absolute','filt_zscored', 'raw_absolute','raw_zscored','time_zscored'};
test2 = isfield(cell_metrics.waveforms,field2remove);
cell_metrics.waveforms = rmfield(cell_metrics.waveforms,field2remove(test2));
end

% cleaning cell_metrics.general.session & cell_metrics.general.animal
field2remove = {'name'};
Expand Down
16 changes: 13 additions & 3 deletions calc_CellMetrics/ce_MonoSynConvClick.m
Original file line number Diff line number Diff line change
Expand Up @@ -147,12 +147,20 @@

% Create CCGs (including autoCG) for all cells
disp('Generating CCGs')
tic
[ccgR1,tR] = CCG(spiketimes,double(spikeIDs(:,3)),'binSize',binSize,'duration',duration,'Fs',1/sr);
toc
ccgR = nan(size(ccgR1,1),nCel,nCel);
ccgR(:,1:size(ccgR1,2),1:size(ccgR1,2)) = ccgR1;


% EAJ added 3/21/2023: also make CCGs binned by time
ccgR_bins = 0:20:spiketimes(end);
ccgR_bins = [ccgR_bins, spiketimes(end)];
for b = 1:length(ccgR_bins)-1
st_binned = spiketimes(spiketimes>ccgR_bins(b) & spiketimes<=ccgR_bins(b+1));
[ccgR1,~] = CCG(st_binned,double(spikeIDs(:,3)),'binSize',binSize,'duration',duration,'Fs',1/sr);
ccgR_temp = nan(size(ccgR1,1),nCel,nCel);
ccgR_temp(:,1:size(ccgR1,2),1:size(ccgR1,2)) = ccgR1;
ccgR_binned{b} = ccgR_temp;
end

% get CI for each CCG
Pval=nan(length(tR),nCel,nCel);
Expand Down Expand Up @@ -294,6 +302,8 @@

% Creating output structure
mono_res.ccgR = ccgR;
mono_res.ccgR_binned = ccgR_binned;
mono_res.ccgR_bins = ccgR_bins;
% mono_res.Pval = Pval;
% mono_res.prob = prob;
% mono_res.prob_noncor = ccgR./permute(repmat(nn2,1,1,size(ccgR,1)),[3 1 2]);
Expand Down
4 changes: 2 additions & 2 deletions calc_CellMetrics/sessionTemplate.m
Original file line number Diff line number Diff line change
Expand Up @@ -230,7 +230,7 @@
% Kilosort
% % % % % % % % % % % % % % % % % % % % % % % % % % % %
if ~exist('relativePath','var')
kiloSortFolder = dir('Kilosort_*');
kiloSortFolder = dir('imec0_ks');
% Extract only those that are directories.
kiloSortFolder = kiloSortFolder(kiloSortFolder.isdir);
if ~isempty(kiloSortFolder)
Expand All @@ -239,7 +239,7 @@
relativePath = ''; % Relative path to the clustered data (here assumed to be the basepath)
end
end
rezFile = dir(fullfile(basepath,relativePath,'rez*.mat'));
rezFile = dir(fullfile(basepath,relativePath,'rez2.mat'));
if ~isempty(rezFile)
rezFile =fullfile(rezFile.folder,rezFile.name);
session = loadKiloSortMetadata(session,rezFile);
Expand Down
39 changes: 39 additions & 0 deletions run_CellExplorer.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,39 @@
%% set session
animal = 'BaggySweatpants';
session_name = '20220214_BaggySweatpants_DY01';
basepath = ['Z:\WT_Sequences\2022_winter\Preprocessed_Data\Spikes\g1\',animal,'\Ecephys\',...
session_name,'\catgt_',session_name,'_g1\',session_name,'_g1_imec0'];
cd(basepath)

%% calculate
tic
% % set params
session = sessionTemplate(basepath);
session.extracellular.fileName = [session_name,'_g1_tcat.imec0.ap.bin'];
session.extracellular.srLfp = 2500;
% session = gui_session(session);
% validateSessionStruct(session);

% generate metrics file
cell_metrics = ProcessCellMetrics('session', session,...
'metrics',{'monoSynaptic_connections'},...
'includeInhibitoryConnections',true,...
'manualAdjustMonoSyn',false,...
'getWaveformsFromDat',false,...
'showWaveforms',false,...
'sessionSummaryFigure',false,...
'showGUI',false);
toc

% %% view
% cell_metrics = CellExplorer('metrics',cell_metrics);
%
% %% view from file
% animal = 'BaggySweatpants';
% session_name = '20220214_BaggySweatpants_DY01';
% basepath = ['Z:\WT_Sequences\2022_winter\Preprocessed_Data\Spikes\g1\',animal,'\Ecephys\',...
% session_name,'\catgt_',session_name,'_g1\',session_name,'_g1_imec0'];
% cd(basepath)
% session = loadSession;
% cell_metrics = loadCellMetrics;
% cell_metrics = CellExplorer('metrics',cell_metrics);
23 changes: 23 additions & 0 deletions run_CellExplorer_Masuda2023.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
addpath(genpath("C:\Users\Niflheim\Documents\GitHub\External\CellExplorer"))
% set session
basepath = "\\oak-smb-giocomo.stanford.edu\groups\giocomo\export\data\Projects\JohnKei_NPH3";
load('\\oak-smb-giocomo.stanford.edu\groups\giocomo\fkmasuda\fkm_analysis\EAJ_revisions\wt_ket_sessions.mat')
for s = 1:length(wt_ket_sessions)
datadir = sprintf('%s\\%s\\%s_keicontrasttrack_ketamine1_g0\\%s_keicontrasttrack_ketamine1_g0_imec0',...
basepath,wt_ket_sessions{s,1},wt_ket_sessions{s,2},wt_ket_sessions{s,2});
cd(datadir)

% set params
session = sessionTemplate(datadir);
session.extracellular.fileName = [wt_ket_sessions{s,2},'_keicontrasttrack_ketamine1_g0_t0.imec0.ap.bin'];

% generate metrics file
cell_metrics = ProcessCellMetrics('session', session,...
'metrics',{'monoSynaptic_connections'},...
'includeInhibitoryConnections',true,...
'manualAdjustMonoSyn',false,...
'getWaveformsFromDat',false,...
'showWaveforms',false,...
'sessionSummaryFigure',false,...
'showGUI',false);
end

0 comments on commit 8da690e

Please sign in to comment.