diff --git a/Analysis/script_loop_through_data_folders.m b/Analysis/script_loop_through_data_folders.m new file mode 100644 index 00000000..a92ce5ab --- /dev/null +++ b/Analysis/script_loop_through_data_folders.m @@ -0,0 +1,67 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_types = {'monopolar', 'bipolar'}; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +t_window = [-2.5, 2.5]; +event_list = {'cueOn', 'centerIn', 'tone', 'centerOut' 'sideIn', 'sideOut', 'foodClick', 'foodRetrieval'}; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.RatID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + for i_lfptype = 1 : length(lfp_types) + + lfp_type = lfp_types{i_lfptype}; + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s', session_name, event_name) + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + + for i_channel = 1 : num_channels + + + end + end + end + end + +end \ No newline at end of file diff --git a/Analysis/script_rename_scalos_to_monopolar.m b/Analysis/script_rename_scalos_to_monopolar.m new file mode 100644 index 00000000..87d1a87a --- /dev/null +++ b/Analysis/script_rename_scalos_to_monopolar.m @@ -0,0 +1,80 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_types = {'monopolar', 'bipolar'}; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +t_window = [-2.5, 2.5]; +event_list = {'cueOn', 'centerIn', 'tone', 'centerOut' 'sideIn', 'sideOut', 'foodClick', 'foodRetrieval'}; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + for i_lfptype = 1 : length(lfp_types) + + lfp_type = lfp_types{i_lfptype}; + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s', session_name, event_name) + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + if ~exist(scalo_folder, 'dir') + continue + end + + for i_channel = 1 : num_channels + [shank_num, site_num] = get_shank_and_site_num(probe_lfp_type, i_channel); + + old_scalo_name = sprintf('%s_scalos_%s_shank%02d_site%02d.mat',session_name, event_name, shank_num, site_num); + old_scalo_name = fullfile(scalo_folder, old_scalo_name); + + if exist(old_scalo_name, 'file') + new_scalo_name = sprintf('%s_scalos_%s_%s_shank%02d_site%02d.mat',session_name, lfp_type, event_name, shank_num, site_num); + new_scalo_name = fullfile(scalo_folder, new_scalo_name); + + movefile(old_scalo_name, new_scalo_name) + end + + end + end + end + end + +end \ No newline at end of file diff --git a/Analysis/sum_session_correct_trials.m b/Analysis/sum_session_correct_trials.m new file mode 100644 index 00000000..d30b3c95 --- /dev/null +++ b/Analysis/sum_session_correct_trials.m @@ -0,0 +1,7 @@ +AllTrials = vertcat(trials.correct); +tone = vertcat(trials.tone); +C = [AllTrials tone]; +indices = find(C(:,1)==0); +C(indices,:) = []; +tone2 = sum( C(:,2)==2 ); +tone1 = sum( C(:,2)==1 ); \ No newline at end of file diff --git a/Analysis/trial_scalograms.m b/Analysis/trial_scalograms.m new file mode 100644 index 00000000..26a96539 --- /dev/null +++ b/Analysis/trial_scalograms.m @@ -0,0 +1,20 @@ +function [event_related_scalos, f, coi] = trial_scalograms(event_triggered_lfps, fb) +%UNTITLED2 Summary of this function goes here +% INPUTS +% event_triggered_lfps - num_trials x samples_per_event array +% Detailed explanation goes here + +num_trials = size(event_triggered_lfps, 1); +samples_per_event = size(event_triggered_lfps, 2); + +f = centerFrequencies(fb); +num_freqs = length(f); + +real_event_related_scalos = zeros(num_trials, num_freqs, samples_per_event); +event_related_scalos = complex(real_event_related_scalos, 0); + +for i_trial = 1 : num_trials + [event_related_scalos(i_trial, :, :), ~, coi] = wt(fb, squeeze(event_triggered_lfps(i_trial, :))); +end + +end \ No newline at end of file diff --git a/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.asv b/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.asv new file mode 100644 index 00000000..c658760b --- /dev/null +++ b/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.asv @@ -0,0 +1,61 @@ +% script create_and_save_trials_struct_Intan + +% loop through all rats with valid choice task data, create trials +% structure for each session, save in the processed data folder to save +% time + +parent_directory = 'Z:\data\ChoiceTask\'; + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_processed_dir = fullfile(session_dirs(i_session).folder, session_name); + cur_rawdata_dir = fullfile(rawdata_folder, session_name); + cd(cur_processed_dir) + + lfp_fname = strcat(session_name, '_lfp.mat'); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_file) + continue + end + + trials_name = sprintf('%s_trials.mat', session_name); + trials_name = fullfile(cur_processed_dir, trials_name); + if exist(trial_name, 'file') + + continue + end + + rawdata_ephys_folder = get_rawdata_ephys_folder(rawdata_folder,session_name); + nexData = extractEventsFromIntanSystem(rawdata_ephys_folder); + + log_file = find_log_file(session_name, parent_directory); + logData = readLogData(log_file); + + sprintf('loaded logData and nexData for %s', session_name) + + trials = createTrialsStruct_simpleChoice_Intan( logData, nexData ); + + + + end + +end \ No newline at end of file diff --git a/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.m b/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.m new file mode 100644 index 00000000..11270a43 --- /dev/null +++ b/ChoiceTask_Intan_behavior_analysis/create_and_save_trials_struct_Intan.m @@ -0,0 +1,83 @@ +% script create_and_save_trials_struct_Intan + +% loop through all rats with valid choice task data, create trials +% structure for each session, save in the processed data folder to save +% time + +parent_directory = 'Z:\data\ChoiceTask\'; + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_processed_dir = fullfile(session_dirs(i_session).folder, session_name); + cur_rawdata_dir = fullfile(rawdata_folder, session_name); + cd(cur_processed_dir) + + lfp_fname = strcat(session_name, '_lfp.mat'); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_file) + continue + end + + trials_name = sprintf('%s_trials.mat', session_name); + trials_name = fullfile(cur_processed_dir, trials_name); + if exist(trials_name, 'file') + % skip if already calculated + sprintf('trials structure already calculated for %s', session_name) + continue + end + + rawdata_ephys_folder = get_rawdata_ephys_folder(rawdata_folder,session_name); + % check that the digitalIn file exists - was missing from some + % early sessions + digitalin_fname = fullfile(rawdata_ephys_folder, 'digitalin.dat'); + if ~exist(digitalin_fname, 'file') + sprintf('no digital input file found for %s', session_name) + continue + end + + nexData = extractEventsFromIntanSystem(rawdata_ephys_folder); + if isempty(nexData) + % something was wrong with the analog/digital input files from + % the intan system + sprintf('nexData could not be generated for %s', session_name) + continue + end + + log_file = find_log_file(session_name, parent_directory); + logData = readLogData(log_file); + + sprintf('loaded logData and nexData for %s', session_name) + + try + trials = createTrialsStruct_simpleChoice_Intan( logData, nexData ); + catch ME + if strcmp(ME.identifier, 'lognexmerge:lognexmismatch') + sprintf('mismatch between log and nex files for %s', session_name) + continue + end + end + + save(trials_name, 'trials'); + + end + +end \ No newline at end of file diff --git a/ChoiceTask_Intan_behavior_analysis/extractEventsFromIntanSystem.m b/ChoiceTask_Intan_behavior_analysis/extractEventsFromIntanSystem.m index acb32328..65a26512 100644 --- a/ChoiceTask_Intan_behavior_analysis/extractEventsFromIntanSystem.m +++ b/ChoiceTask_Intan_behavior_analysis/extractEventsFromIntanSystem.m @@ -28,6 +28,10 @@ digital_word = readIntanDigitalFile(dig_in_file); toc tic +if ~isfield(intan_info, 'board_adc_channels') + nexData = []; + return +end ADC_signals = readIntanAnalogFile(analog_in_file,intan_info.board_adc_channels); toc diff --git a/ChoiceTask_Intan_behavior_analysis/intan2nex.m b/ChoiceTask_Intan_behavior_analysis/intan2nex.m new file mode 100644 index 00000000..cf2f1f34 --- /dev/null +++ b/ChoiceTask_Intan_behavior_analysis/intan2nex.m @@ -0,0 +1,322 @@ +function nexData = intan2nex(dig_in,analog_in,intan_info,varargin) +% +% usage: nexData = intan2nex(dig_in,analog_in,intan_info,varargin) +% +% function to read in timestamps from intan digital and analog in lines and +% convert to a .nex file with events for the choice task +% +% INPUTS: +% dig_in - name of digital input file recorded with the intan system. +% These are generally digital lines from the NI boards running the +% choice task. +% analog_in - name of analog input file recorded with the intan system. +% These are generally digital lines from the NI boards running the +% choice task, plugged into analog lines if we run out of digital +% lines. Note that as of 3/26/2020, the line names are hard-coded in +% the dig_linenames = ... and analog_linenames = ... lines below. +% intan_info - data structure containing info read in from the rhd file +% .XXX - +% TO ADD HERE: DESCRIPTION OF THE INFORMATION IN THE INTAN_INFO FILE, +% WHICH I THINK IS WHAT'S READ IN FROM THE .RHD FILE +% +% VARARGS: +% 'writefile' - whether or not to write the nex formatted data to a new +% file; default filename is the original filename with .nex appended +% - that is, should be "XXXX.box.nex" where XXXX is the original base +% filename. NOT SURE IF WE STILL NEED THIS +% +% OUTPUTS: +% nexData - nex data structure +% nexData.version - file version +% nexData.comment - file comment +% nexData.tbeg - beginning of recording session (in seconds) +% nexData.tend - end of resording session (in seconds) +% +% nexData.neurons - array of neuron structures +% neurons.name - name of a neuron variable +% neurons.timestamps - array of neuron timestamps (in seconds) +% to access timestamps for neuron 2 use {n} notation: +% nexData.neurons{2}.timestamps +% +% nexData.events - array of event structures +% event.name - name of neuron variable +% event.timestamps - array of event timestamps (in seconds) +% to access timestamps for event 3 use {n} notation: +% nexData.events{3}.timestamps +% +% nexData.intervals - array of interval structures +% interval.name - name of neuron variable +% interval.intStarts - array of interval starts (in seconds) +% interval.intEnds - array of interval ends (in seconds) +% +% nexData.waves - array of wave structures +% wave.name - name of neuron variable +% wave.NPointsWave - number of data points in each wave +% wave.WFrequency - A/D frequency for wave data points +% wave.timestamps - array of wave timestamps (in seconds) +% wave.waveforms - matrix of waveforms (in milliVolts), each +% waveform is a vector +% +% nexData.contvars - array of contvar structures +% contvar.name - name of neuron variable +% contvar.ADFrequency - A/D frequency for data points +% +% continuous (a/d) data come in fragments. Each fragment has a timestamp +% and an index of the a/d data points in data array. The timestamp corresponds to +% the time of recording of the first a/d value in this fragment. +% +% contvar.timestamps - array of timestamps (fragments start times in seconds) +% contvar.fragmentStarts - array of start indexes for fragments in contvar.data array +% contvar.data - array of data points (in milliVolts) +% +% nexData.markers - array of marker structures +% marker.name - name of marker variable +% marker.timestamps - array of marker timestamps (in seconds) +% marker.values - array of marker value structures +% marker.value.name - name of marker value +% marker.value.strings - array of marker value strings +% +% +% make it optional to supply a file or list of files; if none supplied, +% then use uigetfile +% +% change log: +% 8/1/12: found that sometimes you get "pops" from high to low to high, +% but can also get them low to high to low. Added code to take +% account of this possibility. +% 03/25/2020 - this may not be necessasry anymore + +analog_thresh = 2; % volts + +writeFile = true; + +% names of what each digital line represents. I think this is correct as of +% 3/25/2020. -DL +dig_linenames = {'', 'foodport', 'food', 'nose5', 'houselight', 'nose4', ... + 'cue5', 'nose3', 'cue4', 'nose2', 'cue3', 'nose1', ... + 'cue2', 'gotrial', 'cue1', 'camframe'}; +% names of what each analog line represents. I think this is correct as of +% 3/25/2020. -DL +analog_linenames = {'tone1','tone2','whitenoise'}; + +% concatenate the linenames so we have a full list +all_linenames = [dig_linenames,analog_linenames]; + +for iarg = 1 : 2 : nargin - 3 + switch lower(varargin{iarg}) + case 'writefile' + writeFile = varargin{iarg + 1}; + case 'analog_thresh' + analog_thresh = varargin{iarg + 1}; + end +end + +num_dig_lines = length(dig_linenames); +% num_analog_lines = length(analog_linenames); + + +% numLines = 16; % 16 digital in lines on the intan system. doesn't count the analog in lines being used as digital lines +numBytesPerSample = 2; % assume uint16's +nSamples = length(dig_in); +totalBytes = nSamples * numBytesPerSample; + +% pull in the sample rate from the intan_info structure +nSampleRate = max(intan_info.frequency_parameters.board_adc_sample_rate,... + intan_info.frequency_parameters.board_dig_in_sample_rate); % these should always be the same +maxNumEvents = 20000; % will choke if there are more than this many events on a line + +% Set up the NEX file data structure +nexData.version = 100; +nexData.comment = 'Converted by intan2nex.m. Dan Leventhal and Jen Magnusson, 2020'; +nexData.freq = nSampleRate; +nexData.tbeg = 0; +nexData.tend = totalBytes/(numBytesPerSample*nSampleRate); +nexData.events = {}; +on.events = {}; +off.events = {}; + +dataOffset = 0; + +usedLines = []; +numUsedLines = 0; +num_usedDigLines = 0; +numLines = length(all_linenames); +for i=1:numLines + if strcmp(all_linenames{i}, ''); continue; end + numUsedLines = numUsedLines + 1; + if i <= num_dig_lines + num_usedDigLines = num_usedDigLines + 1; + end + usedLines(end+1) = i; + + % create nexData.events for each time a line goes "high" or "low". Note + % that transitions from high to low are actually lights turning on, + % nose ports becoming occupied, etc. + + % COMMENTS BELOW NOT RELEVANT RIGHT NOW (3/25/2020) BECAUSE WE'RE NOT USING + % OPTOGENETICS IN THESE EXPERIMENTS... + % this chunk of code will name the events differently for laser + % and shutter lines. 'On' events are generally going from high + % to low (ie, nose-pokes, lighting cue lights, etc) because the + % relays to the MedAssociates boxes are inverted. Laser "on" + % events, however, are indicated by low to high transitions. + % The laser shutter opens on a high value and closes on a low + % value. -DL 02/03/2012 + if contains(lower(all_linenames{i}), 'laser') + % this is a laser line + nexData.events{end+1}.name = [all_linenames{i} 'Off']; + nexData.events{end}.timestamps = {}; + + nexData.events{end+1}.name = [all_linenames{i} 'On']; + nexData.events{end}.timestamps = {}; + + elseif contains(lower(all_linenames{i}), 'video') + % this is a video trigger line + nexData.events{end+1}.name = [all_linenames{i} 'Off']; + nexData.events{end}.timestamps = {}; + + nexData.events{end+1}.name = [all_linenames{i} 'On']; + nexData.events{end}.timestamps = {}; + + elseif contains(lower(all_linenames{i}), 'shutter') + % this is a shutter line + nexData.events{end+1}.name = [all_linenames{i} 'Closed']; + nexData.events{end}.timestamps = {}; + + nexData.events{end+1}.name = [all_linenames{i} 'Open']; + nexData.events{end}.timestamps = {}; + + elseif contains(lower(all_linenames{i}), 'nose') + % this is a nose in/out line + nexData.events{end + 1}.name = [all_linenames{i} 'In']; + nexData.events{end}.timestamps = {}; + + nexData.events{end + 1}.name = [all_linenames{i} 'Out']; + nexData.events{end}.timestamps = {}; + + else + % all other channels + nexData.events{end+1}.name = [all_linenames{i} 'On']; + nexData.events{end}.timestamps = {}; + + nexData.events{end+1}.name = [all_linenames{i} 'Off']; + nexData.events{end}.timestamps = {}; + end + +end + +allLines = false(numUsedLines, nSamples); +allOnEvents = zeros(numUsedLines, maxNumEvents); % maximum 20000 events or so +allOffEvents = zeros(numUsedLines, maxNumEvents); +onEventIdx = ones(numUsedLines,1); +offEventIdx = ones(numUsedLines,1); + +% digital_data = readIntanDigitalFile(dig_in); + +% CREATE A LOOP TO READ IN BLOCKS OF DATA AND CONVERT TO .NEX EVENTS +% Collect line values from every line we care about +for i=1:numUsedLines + if i <= num_usedDigLines + allLines(i,:) = bitget(dig_in, usedLines(i)); % dig_in should really be the array of digital inputs read in + else + analog_idx = i - num_usedDigLines; + allLines(i,:) = (analog_in(analog_idx,:) > analog_thresh); % analog_in should be data read from the analog_in file + end +end + + % Sometimes all the lines are low at the same time + % This is incorrect, and we're going to fix it right here + % in software. + % NOTE: this doesn't explicitly detect "pops", where + % the lines go low then high again, but instead + % just detects whether all lines are low. Should be fine, I _think_... + % above is from AW. + + % this needs to be changed because it's not always true that + % ALL lines "pop" simultaneously. Rewrite the algorithm to + % detect transitions to "low" and then immediately back to + % "high". This happens one SOME but not all lines every time a + % digital port read is initiated in LV. Fucking annoying. There + % may be a way to fix this with a clever combination of + % resistors, but I think this software work-around will be OK. + + % hopefully, this isn't a problem in the Intan system, so I commented it + % out - DL, 12/2019 + +% for iLine = 1 : numLines +% highIdx = find([lastSamp(iLine), allLines(iLine,:)] == 1); +% % compute spacing between high values along each line +% if isempty(highIdx); continue; end; +% +% highDiff = diff(highIdx); +% popIdx = (highDiff == 2); % find indices of events where line went from high to low to high in 3 samples +% +% % 8/1/12 DL%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% lowIdx = find([lastSamp(iLine), allLines(iLine,:)] == 0); +% lowDiff = diff(lowIdx); +% inversePopIdx = (lowDiff == 2); % find indices of events where line went from low to high to low in 3 samples +% % 8/1/12 DL%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% +% allLines(iLine, highIdx(popIdx)) = 1; +% allLines(iLine, lowIdx(inversePopIdx)) = 0; % 8/1/12 DL +% end +% +% % below is left-over from when a line pop was considered when +% % ALL lines go low simultaneously +% % for i=1:length(idx) +% % +% % allLines(:,idx(i)) = allLines(:,idx(i)-1); +% % end +% +% Now we'll extract on and off times +for i=1:numUsedLines + % A channel is considered to have flipped "on" when it goes from 1 to 0 + % A channel is considered to have flipped "off" when it goes from 0 to 1 + + % added lastSamp into the calculation below to account for + % the possibility that a transition takes place exactly + % where a new block of values is loaded. -DL 20110627 + onTimes = find( diff(allLines(i,:)) == -1 )/nSampleRate; + offTimes = find( diff(allLines(i,:)) == 1 )/nSampleRate; + numOn = length(onTimes); + numOff = length(offTimes); + allOnEvents(i,onEventIdx(i):onEventIdx(i)+numOn-1) = onTimes; + allOffEvents(i,offEventIdx(i):offEventIdx(i)+numOff-1) = offTimes; + onEventIdx(i) = onEventIdx(i) + numOn; + offEventIdx(i) = offEventIdx(i) + numOff; +end + +% lastSamp = allLines(:, end); +% + + +% Stick in the events +for i=1:numUsedLines + nexData.events{i*2-1}.timestamps = unique(allOnEvents(i,1:onEventIdx(i)-1)'); + nexData.events{i*2}.timestamps = unique(allOffEvents(i,1:offEventIdx(i)-1)'); +end +nexData.events = nexData.events'; + +% Make intervals +nexData.intervals = {}; +for i=1:numUsedLines + + intStarts = nexData.events{i*2-1}.timestamps; + intEnds = nexData.events{i*2}.timestamps; + + if isempty(intEnds) || isempty(intStarts); continue; end + + intEnds = intEnds( intEnds > intStarts(1) ); + + if isempty(intEnds); continue; end + intStarts = intStarts( intStarts < intEnds(end) ); + + nexData.intervals{end+1}.name = all_linenames{usedLines(i)}; + nexData.intervals{end}.intStarts = intStarts; % notice we used end+1 above. that created the entry, and now we can just refer to it as "end" + nexData.intervals{end}.intEnds = intEnds; +end +nexData.intervals = nexData.intervals'; % stupid, STUPID requirement for writeNexFile... + +% if writeFile +% writeNexFile(nexData, [filename '.nex']); +% end \ No newline at end of file diff --git a/Helpers/display_scalogram.m b/Helpers/display_scalogram.m new file mode 100644 index 00000000..e7780024 --- /dev/null +++ b/Helpers/display_scalogram.m @@ -0,0 +1,48 @@ +function [outputArg1,outputArg2] = display_scalogram(scalo, t_window, fb, varargin) +%UNTITLED Summary of this function goes here +% Detailed explanation goes here + +cmap = 'jet'; +h_ax = 0; +c_lim = [0, 0]; + +f_ticks = [1, 10, 20, 50, 100]; + +for i_arg = 4 : 2 : nargin - 1 + switch lower(varargin{iarg}) + case 'colormap' + cmap = varargin{iarg + 1}; + case 'ax' + h_ax = varargin{iarg + 1}; + case 'fticks' + f_ticks = varargin{iarg + 1}; + case 'clim' + c_lim = varargin{iarg + 1}; + end +end + +if h_ax == 0 + figure; + h_ax = gca; +end + +num_samples = size(scalo, 2); + +t = linspace(t_window(1), t_window(2), num_samples); +f = centerFrequencies(fb); + +axes(h_ax) + +surface(t, f, scalo); + +axis tight +shading flat + +if any(c_lim ~= 0) + set(gca,'clim',c_lim) +end +set(gca,'ytick',f_ticks); +set(gca,'yscale','log') + +colormap(cmap) +end \ No newline at end of file diff --git a/Helpers/extract_perievent_data.m b/Helpers/extract_perievent_data.m new file mode 100644 index 00000000..db029d00 --- /dev/null +++ b/Helpers/extract_perievent_data.m @@ -0,0 +1,29 @@ +function perievent_data = extract_perievent_data(ephys_data, trials, event_name, t_window, Fs) +%UNTITLED2 Summary of this function goes here +% INPUTS +% ephys_data - num_channels x num_samples array +% ts - vector of timestamps with times in seconds +% +% Detailed explanation goes here + +ts = ts_from_trials(trials, event_name); +num_trials = length(ts); +num_channels = size(ephys_data, 1); +total_samples = size(ephys_data, 2); + +samp_window = round(t_window * Fs); +center_samps = round(ts * Fs); +samp_windows = center_samps + samp_window; +samps_per_window = range(samp_window) + 1; + +perievent_data = zeros(num_trials, num_channels, samps_per_window); + +for i_trial = 1 : num_trials + + if samp_windows(i_trial, 1) < 1 || samp_windows(i_trial, 2) > total_samples + % if window starts before start of recording or ends after end of + % recording, skip + continue + end + perievent_data(i_trial, :, :) = ephys_data(:, samp_windows(i_trial, 1) : samp_windows(i_trial, 2)); +end \ No newline at end of file diff --git a/Helpers/extract_trials_by_features.m b/Helpers/extract_trials_by_features.m new file mode 100644 index 00000000..602fde58 --- /dev/null +++ b/Helpers/extract_trials_by_features.m @@ -0,0 +1,63 @@ +function [valid_trials, valid_trial_flags] = extract_trials_by_features(trials, trialfeatures) +% INPUTS +% trials +% trialfeatures - string containing trial features to extract. If any of +% the following strings are containes in 'trialfeatures', the +% following will be extracted. Can do this in any combination +% 'correct' - extracts correct trials +% 'wrong' - extracts incorrect trials +% 'moveright' - extracts trials in which rat moved right +% 'moveleft' - extracts trials in which rat moved left +% 'cuedleft' - extracts trials in which tone prompted rat to move +% left +% 'cuedright' - extracts trials in which tone prompted rat to move +% right +% 'falsestart' - extracts false start trials +% Detailed explanation goes here + +num_trials = length(trials); +valid_trial_flags = true(num_trials, 1); +if contains(lower(trialfeatures), 'correct') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'correct', true); +end + +if contains(lower(trialfeatures), 'wrong') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'correct', false); +end + +if contains(lower(trialfeatures), 'moveright') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'movementDirection', 1); +end + +if contains(lower(trialfeatures), 'moveleft') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'movementDirection', 2); +end + +if contains(lower(trialfeatures), 'cuedleft') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'tone', 1); +end + +if contains(lower(trialfeatures), 'cuedright') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'tone', 2); +end + +if contains(lower(trialfeatures), 'falsestart') + valid_trial_flags = valid_trial_flags & find_trials_by_field(trials, 'falseStart', 1); +end + +valid_trials = trials(valid_trial_flags); + +end + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +function is_valid_trial = find_trials_by_field(trials, field, value) + +num_trials = length(trials); +is_valid_trial = false(num_trials, 1); + +for i_trial = 1 : num_trials + is_valid_trial(i_trial) = (trials(i_trial).(field) == value); +end + +end \ No newline at end of file diff --git a/Helpers/extract_trials_by_outcome.m b/Helpers/extract_trials_by_outcome.m new file mode 100644 index 00000000..b1237f09 --- /dev/null +++ b/Helpers/extract_trials_by_outcome.m @@ -0,0 +1,28 @@ +function [outputArg1,outputArg2] = extract_trials_by_outcome(trials, trialoutcome) +% INPUTS +% trials +% trialoutcome +% Detailed explanation goes here + +valid_trial_idx = +if contains(lower(trialoutcome), 'correct') + +end + +if contains(lower(trialoutcome), 'moveright') + +end + +if contains(lower(trialoutcome), 'moveleft') + +end +outputArg1 = inputArg1; +outputArg2 = inputArg2; +end + + +function is_valid_trial = find_trials_by_field(trials, field, value) + +num_trials = length(trials); + +for i_trial = 1 : num_trials diff --git a/Helpers/find_physiology_data.m b/Helpers/find_physiology_data.m new file mode 100644 index 00000000..663382ca --- /dev/null +++ b/Helpers/find_physiology_data.m @@ -0,0 +1,43 @@ +function physiology_folder = find_physiology_data(session_dir) + +physiology_folder = ''; + +[~, session_name, ~] = fileparts(session_dir); +str_parts = split(session_name, '_'); +ratID = str_parts{1}; + + +if length(str_parts) == 1 + % workaround for a poorly named folder + return +end +if strcmpi('str_parts{2}', 'check') + % workaround for a poorly named pair of folders + return +end + +session_datestr = str_parts{2}(1:8); % assume yyyymmdd date format + +rat_date = strcat(ratID, '_', session_datestr); + +subdirs = dir(fullfile(session_dir, strcat(rat_date, '*'))); + +num_subdirs = length(subdirs); + +for i_dir = 1 : num_subdirs + + cur_path = fullfile(session_dir, subdirs(i_dir).name); + if isfolder(cur_path) + + % look for an amplifier.dat and info.rhd file + amp_file = fullfile(cur_path, 'amplifier.dat'); + info_file = fullfile(cur_path, 'info.rhd'); + + if isfile(amp_file) && isfile(info_file) + physiology_folder = cur_path; + end + end + +end + +end \ No newline at end of file diff --git a/Helpers/get_rat_list.m b/Helpers/get_rat_list.m new file mode 100644 index 00000000..d011c39d --- /dev/null +++ b/Helpers/get_rat_list.m @@ -0,0 +1,19 @@ +function [rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list() +% list of rats to be analyzed from thalamic recordings +% update rat_nums (and hopefully not ratnums_with_bad_histo) to indicate +% which rats should be analyzed + +rat_nums = [326, 327, 372, 374, 376, 378, 379, 394, 395, 396, 411, 412, 413, 419, 420, 425, 427, 456, 459, 460, 462, 464, 463, 466, 465, 467, 479]; + +ratnums_with_bad_histo = [374, 376, 396, 413]; + +rats_with_good_histo = 0; +for i_rat = 1 : length(rat_nums) + ratIDs{i_rat} = sprintf('R%04d', rat_nums(i_rat)); + + if ~ismember(rat_nums(i_rat), ratnums_with_bad_histo) + rats_with_good_histo = rats_with_good_histo + 1; + + ratIDs_goodhisto{rats_with_good_histo} = ratIDs{i_rat}; + end +end \ No newline at end of file diff --git a/Helpers/get_shank_and_site_num.m b/Helpers/get_shank_and_site_num.m new file mode 100644 index 00000000..898a4b1c --- /dev/null +++ b/Helpers/get_shank_and_site_num.m @@ -0,0 +1,47 @@ +function [column_num,site_num] = get_shank_and_site_num(probe_type, row_idx) +%UNTITLED7 Summary of this function goes here +% INPUTS +% row_idx = the index of the row in the ordered_lfp array. It should be +% ordered such that we start at shank 1, most dorsal site, then go +% down the shank to site 8 (most ventral sitefor an NN8x8), then +% shank 2 running from dorsal to ventral, etc. For a Cambridge probe, +% it should start at the most dorsal site and go down 16 sites, with +% 4 columns of electrodes total +% +% OUTPUTS +% column_num - number of recording column - could be a single shank, or +% there could be multiple columns of recording sites per shank +% site_num - site going from dorsal to ventral along the column (note, +% might not be a physical site for bipolar). e.g., site 1 is the most +% dorsal, site 2 ventral to site 1, etc, then start over in the next +% column + +probe_parts = split(probe_type, '_'); + +% which style of probe was used? +switch lower(probe_parts{1}) + case 'nn8x8' + % 8 shanks with 8 sites each + sites_per_column = 8; + + case 'assy156' + sites_per_column = 16; + + case 'assy236' + sites_per_column = 16; + +end + +% monopolar vs bipolar +switch lower(probe_parts{2}) + + case 'monopolar' + signals_per_column = sites_per_column; + + case 'bipolar' + signals_per_column = sites_per_column - 1; + +end + +column_num = ceil(row_idx / signals_per_column); +site_num = row_idx - (column_num-1) * signals_per_column; \ No newline at end of file diff --git a/Helpers/plotSpikeRaster_color.m b/Helpers/plotSpikeRaster_color.m index b493c905..ddb22692 100644 --- a/Helpers/plotSpikeRaster_color.m +++ b/Helpers/plotSpikeRaster_color.m @@ -1,22 +1,43 @@ -function plotSpikeRaster_color(xPoints,yPoints,groups,colors,figPos,markerSize) +function plotSpikeRaster_color(unitBehavior,event_idx,groups,colors,figPos,markerSize) +% INPUTS +% xPoints (unitBehavior) - vector of x coordinates at which to put a raster mark +% (timestamps relative to a behavioral event for an entire session) +% yPoints (event_idx) - vector of y coordinates at which to put a raster mark +% xPoints and yPoints should have the same length, and xPoints(1) +% goes with yPoints(1), xPoints(2)-->yPoints(2), etc. +% the unique values of yPoints indicates a single trial, so all +% indicates of yPoints with the same value can be used to find +% xPoints that go with them, and plot them all the same color along +% the same row +% yPoints should be integers representing trial #'s (i.e., 1, 2, 3, +% ...) +% groups - vector containing num_trials (number of unique yPoints) list +% integers, with each integer representing a different group (for +% example, maybe 1 - leftward movement, 2-rightward movement, then +% they will show up as different colors in the raster +% colors - vector of colors corresponding to groups (i.e., colors(1) is +% the color that all members of group 1 will be plotted) +% figPos - figure position +% markerSize - size of the marker to use for each raster point + % % if numel(xPoints) ~= numel(yPoints) || numel(unique(yPoints)) ~= numel(groups) ||... % % numel(unique(groups)) ~= size(colors,1) % % error('check input dimensions'); % % end -yPoints_unique = unique(yPoints); +yPoints_unique = unique(event_idx); if ~isempty(figPos) figure('position',figPos); end -for iTrial = 1:numel(unique(yPoints)) +for iTrial = 1:numel(unique(event_idx)) cur_y = yPoints_unique(iTrial); - pointIdxs = find(yPoints == cur_y); - plot(xPoints(pointIdxs),yPoints(pointIdxs),'.','color',colors(groups(iTrial) == unique(groups),:),'markerSize',markerSize); + pointIdxs = find(event_idx == cur_y); + plot(unitBehavior(pointIdxs),event_idx(pointIdxs),'.','color',colors(groups(iTrial) == unique(groups),:),'markerSize',markerSize); hold on; end -ylim([1 numel(unique(yPoints))]); +ylim([1 numel(unique(event_idx))]); set(gca,'ydir','reverse'); if ~isempty(figPos) xlabel('time (s)'); diff --git a/Helpers/readLogData.m b/Helpers/readLogData.m new file mode 100644 index 00000000..1bf0b885 --- /dev/null +++ b/Helpers/readLogData.m @@ -0,0 +1,98 @@ +function logData = readLogData(fname) +% +% update 08/22/2014 to allow varargin for filename +% +% usage: logData = readLogData( fname ) +% +% INPUTS: +% fname - name of the .log file +% +% OUTPUTS: +% logData - structure with the following fields: +% fileVersion (uint16) - version of the logData file written by LabView +% taskID (uint8) - this is a task identifier (ie, stop-signal, +% go-nogo, simple choice, RL, food-water, etc.). This number +% should match with the code used in the sql database +% taskVersion - version number of the specific task. This differs +% from "fileVersion" in that fileVersion refers to the log file +% structure, while "taskVersion" refers to the version of the +% actual data written for the specific task +% subject - string containing the subject name. Max 10 characters in +% fileVersion 1 +% date - string containing the data formatted as yyyymmdd +% startTime - string containing 24-hour start time based on the +% behavior computer clock, in 'HH:MM' format. +% comment - string containing a session comment (max 1024 characters) +% +% After that is a collection of "header" fields. The names of these +% fields are written into the .log file. These are parameters +% that do not change over the course of a session and are +% task-specific. +% +% After that comes the actual data. The structure data field names +% are also written into the .log file and are task-specific. + +%fname = varargin{1}; //there's only one input, so no need for varargin +bitOrder = 'b'; + +commentLength = 200; +%open the file that was input to read from +fid = fopen(fname, 'r'); + +logData.fileVersion = fread(fid, 1, 'uint16', 0, bitOrder); +if logData.fileVersion > 1000 % indicates that this is not the standard log file for Leventhal lab, could be a legacy from Berke lab + logData = []; + return; +end + +%fill in fields of logData from the file +logData.taskID = fread(fid, 1, 'uint8', 0, bitOrder); +logData.taskVersion = fread(fid, 1, 'uint8', 0, bitOrder); +logData.subject = deblank(fread(fid, 10, '*char')'); +logData.date = fread(fid, 8, '*char')'; +logData.startTime = fread(fid, 5, '*char')'; + +%change positions within the file +fseek(fid, 2 * 1024, 'bof'); + +%get rid of null characters within the comment +logData.comment = deblank(fread(fid, 1024, '*char')'); + +fseek(fid, 3 * 1024, 'bof'); +% read in the header fields +fullHeaderString = deblank(fread(fid, 1024, '*char')'); +dlmIdx = findstr(fullHeaderString, ','); % find the location of all the commas +dlmIdx = [0, dlmIdx, length(fullHeaderString)+1]; +numHeaderFields = length(dlmIdx) - 1; +headerFieldNames = cell(1, numHeaderFields); + +for iField = 1 : numHeaderFields + headerFieldNames{iField} = fullHeaderString(dlmIdx(iField)+1:dlmIdx(iField+1)-1); +end + +fseek(fid, 4 * 1024, 'bof'); +% read in the values for the header fields +for iField = 1 : numHeaderFields + logData.(headerFieldNames{iField}) = fread(fid, 1, 'double', 0, bitOrder); +end + +fseek(fid, 5 * 1024, 'bof'); +% read in the data fields +fullDataString = deblank(fread(fid, 1024, '*char')'); +dlmIdx = findstr(fullDataString, ','); % find the location of all the commas +dlmIdx = [0, dlmIdx, length(fullDataString)+1]; +numDataFields = length(dlmIdx) - 1; +dataFieldNames = cell(1, numDataFields); +for iField = 1 : numDataFields + dataFieldNames{iField} = fullDataString(dlmIdx(iField)+1:dlmIdx(iField+1)-1); +end + +fseek(fid, 6 * 1024, 'bof'); +% read in the actual data +data = fread(fid, [numDataFields, inf], 'double', 0, bitOrder); +if ~isempty(data) + for iField = 1 : numDataFields + logData.(dataFieldNames{iField}) = data(iField, :)'; + end +end +fclose(fid); \ No newline at end of file diff --git a/Helpers/read_Jen_xls_summary.m b/Helpers/read_Jen_xls_summary.m new file mode 100644 index 00000000..7ff7e2e5 --- /dev/null +++ b/Helpers/read_Jen_xls_summary.m @@ -0,0 +1,16 @@ +function xls_data = read_Jen_xls_summary(filename, sheetname) +%UNTITLED8 Summary of this function goes here +% Detailed explanation goes here + +% update this function when need information from the other spreadsheets + +switch sheetname + case 'probe_type' + xls_data = readtable(filename,'filetype','spreadsheet',... + 'sheet','probe_type',... + 'texttype','string'); + otherwise + xls_data = ''; +end + +end \ No newline at end of file diff --git a/Helpers/ts_from_trials.m b/Helpers/ts_from_trials.m new file mode 100644 index 00000000..3f822976 --- /dev/null +++ b/Helpers/ts_from_trials.m @@ -0,0 +1,17 @@ +function ts = ts_from_trials(trials, event_name) +%UNTITLED6 Summary of this function goes here +% Detailed explanation goes here + +num_trials = length(trials); + +ts = NaN(num_trials, 1); + +for i_trial = 1 : num_trials + + if isfield(trials(i_trial).timestamps, event_name) + ts(i_trial) = trials(i_trial).timestamps.(event_name); + end + +end + +end diff --git a/LFPs/Analysis_code_JM/LFP_Analysis_Workflow.docx b/LFPs/Analysis_code_JM/LFP_Analysis_Workflow.docx index 6906db72..dfc35df6 100644 Binary files a/LFPs/Analysis_code_JM/LFP_Analysis_Workflow.docx and b/LFPs/Analysis_code_JM/LFP_Analysis_Workflow.docx differ diff --git a/LFPs/Analysis_code_JM/PETHS/ATH_JM_peths.m b/LFPs/Analysis_code_JM/PETHS/ATH_JM_peths.m new file mode 100644 index 00000000..260c7e02 --- /dev/null +++ b/LFPs/Analysis_code_JM/PETHS/ATH_JM_peths.m @@ -0,0 +1,66 @@ +function [unitBehavior, event_idx, pethEntries] = ATH_JM_peths(unit_struct, trial_ts) +% creating Alex's as a function + +% INPUTS +% unit_struct - matrix for each neuron, could be for each channel +% trial_ts - trial_ts is an m x n array where m is num_trials and n is +% samples. This is generated looking for a specific time of the trial +% structure (e.g. correct trials). Units should be seconds +% +% OUTPUTS +% unitBehavior - vector containing peri-event timestamps +% event_idx = cell array of vectors containing the index of the event +% around which unitBehavior is calculated (for example, if +% unitBehavior is [-0.5, 0.1, 0.5, -0.6] and event_idx is [1, 1, 1, +% 2], that would mean the first 3 timestamps are around the first +% event in the recording and the last timestamp (-0.6) is from the +% second event. unitBehavior and event_idx can be input to +% plotSpikeRaster_color as xPoints and yPoints, respectively + + +pre=3; +post=3; +binW=0.050; +pethEntries=[]; +% unitBehavior = cell(num_units, 1); % Jen you added this in +byUnit_Spikes4raster=struct('trial',[],'spike_times',[]); +num_units = length(unit_struct); +event_idx = cell(num_units, 1); +for iUnit=1:length(unit_struct) % matrix for each neuron, could be for each channel + unit_spike_times = unit_struct(iUnit).ts;% spike times for each unit + unitBehavior=[]; % matrix of this unit's behavior-related activity + %unitBehavior=nan(length(-pre:binW:post),length(trial_ts)); + %yPoints=[]; + i_event=1; + byUnit_Spikes4raster(iUnit).trial=[]; + byUnit_Spikes4raster(iUnit).spike_times=[]; + for bt=trial_ts %for each behavior instance; I think trial_ts is an m x n array where m is num_trials and n is samples. Is this relative to trial start or each event? + ts_to_append = unit_spike_times(unit_spike_times>(bt-pre) & unit_spike_times<(bt+post))-bt; + % thisN=[thisN, nt(nt>bt-pre & nt(bt-pre) & unit_spike_times<(bt+post))-bt; + % thisN=[thisN, nt(nt>bt-pre & nt(bt-pre) & unit_spike_times<(bt+post))-bt; + + % make sure ts_to_append is a row vector + if iscolumn(ts_to_append) + ts_to_append = ts_to_append'; + end + + % thisN=[thisN, nt(nt>bt-pre & nt trial_ts(iEvent) - pethHalfWidth & sClust < trial_ts(iEvent) + pethHalfWidth); + % subtract starting value to center ts entries + if ~isempty(tsFitCriteria) + pethEntries = [pethEntries sClust(tsFitCriteria) - trial_ts(iEvent)]; + end + end + allTs{totalUnits} = {sClust,unitName}; + totalUnits = totalUnits + 1; + + figure; + histogram(tsFitCriteria,binSize); + xlabel('time (s)'); + ylabel('spikes'); + title([strrep(clusternumbers(iUnit),'_','-')]); + disp(length(pethEntries)) +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/PETHS/utils/create_cluster_timestamps_data_folder.m b/LFPs/Analysis_code_JM/PETHS/utils/create_cluster_timestamps_data_folder.m new file mode 100644 index 00000000..74ea036e --- /dev/null +++ b/LFPs/Analysis_code_JM/PETHS/utils/create_cluster_timestamps_data_folder.m @@ -0,0 +1,10 @@ +function session_cluster_timestamps = create_cluster_timestamps_data_folder(rd_metadata, parent_directory) + +session_clust_timestamps_folder= strcat(rd_metadata.ratID, '-cluster-timestamps'); % cluster_timestamps_data_folder +session_cluster_timestamps = fullfile(parent_directory, rd_metadata.ratID, session_clust_timestamps_folder, rd_metadata.session_name); + +if ~isfolder(session_cluster_timestamps) + + mkdir(session_cluster_timestamps) + +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/extract_LFP_around_timestamps.m b/LFPs/Analysis_code_JM/extract_LFP_around_timestamps.m index bf3b0003..c9975764 100644 --- a/LFPs/Analysis_code_JM/extract_LFP_around_timestamps.m +++ b/LFPs/Analysis_code_JM/extract_LFP_around_timestamps.m @@ -37,11 +37,10 @@ else start_sample = floor(trial_ts(i_ts) * Fs); end_sample = start_sample + samples_per_event - 1; - try current_lfp = ordered_lfp(:, start_sample:end_sample); % changed LFP to ordered_lfp for the calculate_cwt_3D_matrix script catch - keyboard; + keyboard end end diff --git a/LFPs/Analysis_code_JM/extract_trial_ts.m b/LFPs/Analysis_code_JM/extract_trial_ts.m index 99afc6de..0b11b806 100644 --- a/LFPs/Analysis_code_JM/extract_trial_ts.m +++ b/LFPs/Analysis_code_JM/extract_trial_ts.m @@ -1,36 +1,35 @@ -function trial_ts = extract_trial_ts(trials, eventFieldnames) -% +function trial_ts = extract_trial_ts(trials, eventFieldname) % % INPUTS % trials - trials structure - -% % tWindow - two-element array containing length of time (in seconds) -% % before and after each event to extract. Note that both elements of -% % tWindow will be added to the reference event (so usually the first -% % number will be negative, second will be positive) - -% WHY DOES THIS FILE NEED TWINDOW?? DOES NOT MAKE SENSE. - -% eventFieldnames - event fields for which to extract time -% windows +% eventFieldname - name of the event we're interested in % For a correct trial, eventFieldnames are as follows: % cueOn, centerIn, tone, centerOut, sideIn, sideOut, foodClick, foodRetrieval % - to find eventFieldnames for non-correct trials, check the % timestamps column in the trials structure. % % OUTPUTS -% trialRanges - m x n x 2 array, where m is the number of fields being -% analyzed, n is the number of trials, and the last two elements -% contain the start and end time for each window +% trial_ts - vector of timestamps for the given event in this trials +% structure -trial_ts = NaN(numel(eventFieldnames),numel(trials)); -for iField = 1:numel(eventFieldnames) - for iTrial = 1:numel(trials) - try - ts = getfield(trials(iTrial).timestamps,eventFieldnames{iField}); - trial_ts(iField, iTrial) = ts; - catch - % do nothing, filled with NaN - end +trial_ts = NaN(1,numel(trials)); +for iTrial = 1:numel(trials) + try + ts = getfield(trials(iTrial).timestamps, eventFieldname); + trial_ts(iTrial) = ts; + catch end -end \ No newline at end of file +end + + +% trial_ts = NaN(numel(eventFieldnames),numel(trials)); +% for iField = 1:numel(eventFieldnames) +% for iTrial = 1:numel(trials) +% try +% ts = getfield(trials(iTrial).timestamps,eventFieldnames{iField}); +% trial_ts(iField, iTrial) = ts; +% catch +% % do nothing, filled with NaN +% end +% end +% end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data.m b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data.m index 72e1f049..a0269ef9 100644 --- a/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data.m +++ b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data.m @@ -32,9 +32,9 @@ % lists for ratID probe_type NN8x8 = ["R0326", "R0327", "R0372", "R0379", "R0374", "R0376", "R0378", "R0394", "R0395", "R0396", "R0412", "R0413"]; % Specify list of ratID associated with each probe_type ASSY156 = ["R0411", "R0419"]; -ASSY236 = ["R0420", "R0425", "R0427", "R0457"]; +ASSY236 = ["R0420", "R0425", "R0427", "R0457", "R0456", "R0459", "R0460", "R0463", "R0465", "R0466", "R0467", "R0477"]; -sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a', 'R0427_20220920a'}; +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a', 'R0427_20220920a', 'R0427_20220919a'}; sessions_to_ignore1 = {'R0425_20220728_ChVE_220728_112601', 'R0427_20220920_Testing_220920_150255'}; sessions_to_ignore2 = {'R0427_20220908a', 'R0427_20220909a', 'R0427_20220912a','R0427_20220913a', 'R0427_20220914a', 'R0427_20220915a', 'R0427_20220916a'}; % Trying this as a workaround. Code wouldn't skip these two trials. R0425 - 15 hour session and R0427 no data (files didn't save correctly)? @@ -72,10 +72,7 @@ continue; end -% if contains(ratID, NN8x8) || contains(ratID, ASSY156) % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes -% continue; -% end - + % if contains(ratID, 'R0326') || contains(ratID, 'R0372') || contains(ratID, 'R0376')... % || contains(ratID, 'R0374') || contains(ratID, 'R0378') || contains(ratID, 'R0379') % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes % continue; @@ -227,7 +224,7 @@ % Order the lfps here [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); - % Orders the lfps by probe site mapping (double check the single file to remove catches for loading in single data) + % Orders the lfps by probe site mapping (double check the single file to remove catches for loading in single data) %find the column header that matches session_name (session_name and %the column headers should match for each session to pull in the @@ -240,11 +237,10 @@ % write code here to get event_triggered_lfps_ordered % for each separate event -% % create a filename to save the trials_structure + % create a filename to save the trials_structure trials_fname_to_save = char(strcat(session_name, '_', eventFieldnames{i_event}, '_', trialType, '_', 'trials', '.mat')); % add in ''_', sprintf('trial%u.pdf', trial_idx)' should you want to save files individually trials_full_name = fullfile(session_trials_folder, trials_fname_to_save); - % if exist(trials_full_name, 'file') % continue % end diff --git a/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots.m b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots.m index 14a378ba..cab47252 100644 --- a/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots.m +++ b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots.m @@ -4,6 +4,7 @@ intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; rats_with_intan_sessions = find_rawdata_folders(intan_parent_directory); +valid_trials_original_folder = find_trials_struct_original_folders(intan_parent_directory); %% fname = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary\Rat_Information_channels_to_discard.xlsx'; % for channels to ignore based on visualizing in Neuroscope etc @@ -72,10 +73,7 @@ continue; end -% if contains(ratID, NN8x8) || contains(ratID, ASSY156) % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes -% continue; -% end - + % if contains(ratID, 'R0326') || contains(ratID, 'R0372') || contains(ratID, 'R0376')... % || contains(ratID, 'R0374') || contains(ratID, 'R0378') || contains(ratID, 'R0379') % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes % continue; @@ -227,7 +225,7 @@ % Order the lfps here [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); - % Orders the lfps by probe site mapping (double check the single file to remove catches for loading in single data) + % Orders the lfps by probe site mapping (double check the single file to remove catches for loading in single data) %find the column header that matches session_name (session_name and %the column headers should match for each session to pull in the @@ -240,11 +238,10 @@ % write code here to get event_triggered_lfps_ordered % for each separate event -% % create a filename to save the trials_structure + % create a filename to save the trials_structure trials_fname_to_save = char(strcat(session_name, '_', eventFieldnames{i_event}, '_', trialType, '_', 'trials', '.mat')); % add in ''_', sprintf('trial%u.pdf', trial_idx)' should you want to save files individually trials_full_name = fullfile(session_trials_folder, trials_fname_to_save); - % if exist(trials_full_name, 'file') % continue % end @@ -280,110 +277,110 @@ % num_trials = size(event_triggered_lfps, 1); % % catch for num_trials_to_plot < num_trials (this is specifically for R0420_20220714 since it did only 3 trials so can't graph 5 random trials). - if num_trials < num_trials_to_plot - continue; + if num_trials < num_trials_to_plot + continue; + end + + % select trials for plotting at random + trial_idx_to_plot = randperm(num_trials, num_trials_to_plot); % create a catch if num_trials_to_plot < num_trials, continue. + + for i_trial = 1:num_trials_to_plot + trial_idx = trial_idx_to_plot(i_trial); + channel_lfps = squeeze(event_triggered_lfps(trial_idx,:, :)); + % trials, create a 64channel graph for each trial? + + % Plot the data + figure; + num_rows = size(event_triggered_lfps,2); + LFPs_per_shank = num_rows/ 8; % will be 8 for 64 channels, 7 for 56 channels (diff) + for i_row = 1 : num_rows + + y_lim = [-2500, 2500]; + plot_col = ceil(i_row / LFPs_per_shank); + plot_row = i_row - LFPs_per_shank * (plot_col-1); + plot_num = (plot_row-1) * 8 + plot_col; + + subplot(LFPs_per_shank,8,plot_num); + + % sue the trials structure data + % (trials_validchannels_marked.is_channel_valid_ordered) + % to color the actual plots based on good or + % bad sites above the specified threshold (1500mV, check the script to verify) + % switch trials_validchannels_marked(i_trial).is_channel_valid_ordered(i_row) + switch trials_validchannels_marked(trial_idx).is_channel_valid_ordered(i_row) + % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps + case 0 + plot_color = 'r'; % marks bad channels within specified trial + case 1 + plot_color = 'k'; % marks good channels within specified trial + otherwise + plot_color = 'b'; % catch in case the data was not input into the structure + end + + plot_channel_lfps = plot(t, channel_lfps(i_row, :), plot_color); % change to log10 -- plot(f, 10*log10(power_lfps(:,1))) + set(gca, 'ylim',y_lim); + grid on + + % give titles to each subplot + if contains(ratID, NN8x8) % if the ratID is in the list, it'll assign it the correct probe_type for ordering the LFP data correctly + caption = sprintf('NN8x8 #%d', site_order(i_row)); + elseif contains(ratID, ASSY156) + caption = sprintf('ASSY156 #%d', site_order(i_row)); + elseif contains(ratID, ASSY236) + caption = sprintf('ASSY236 #%d', site_order(i_row)); + end + title(caption, 'FontSize', 8); + + if plot_row < LFPs_per_shank + set(gca,'xticklabels',[]) + end + + if plot_col > 1 + set(gca,'yticklabels',[]) + end + + ax = gca; + % This section is coded to color the axes of + % the plots when checking the amplifier.dat + % files 'by eye' using Neuroscope + switch valid_sites_reordered(i_row) % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps + case 0 + ax.XColor = 'r'; % Red % marks bad channels within specified trial + ax.YColor = 'r'; % Red + % ax.ylabel = 'k'; + case 1 + ax.XColor = 'k'; % black % marks good channels within specified trial + ax.YColor = 'k'; % black + % ax.ylabel = 'k'; + case 2 + ax.XColor = 'b'; % blue % marks channels as 'variable' and could be good for portions of the whole amplifier.dat file but bad for others. Thus some channels may be good for only some trials, not all. + ax.YColor = 'b'; + % ax.ylabel = 'k'; + otherwise + ax.XColor = 'b'; % blue % catch in case the data was not input into the structure + ax.YColor = 'b'; + % ax.ylabel = 'k'; + end + end + set(gcf, 'WindowState', 'maximize'); % maximizes the window so that it exports the graphics with appropriate font size + + % Gives the whole plot page a name + A=cell(1,5); + A{1} = ['Subject: ' ratID]; + A{2} = ['Session: ' session_name]; + A{3} = ['Task Level: ' choiceRTdifficulty{logData.taskLevel+1}]; + A{4} = ['Trial Number: ' num2str(trial_idx)]; + A{5} = ['EventFieldname and Trial Type: ' eventFieldnames{i_event} '_' trialType]; + + sgtitle(A, 'Interpreter','none'); + + if i_trial == 1 + exportgraphics(gcf, trials_plot_full_name); + else + exportgraphics(gcf, trials_plot_full_name, 'append', true); %ADD APPEND FLAG HERE + end + close; + end end - - % select trials for plotting at random - trial_idx_to_plot = randperm(num_trials, num_trials_to_plot); % create a catch if num_trials_to_plot < num_trials, continue. - - for i_trial = 1:num_trials_to_plot - trial_idx = trial_idx_to_plot(i_trial); - channel_lfps = squeeze(event_triggered_lfps(trial_idx,:, :)); - % trials, create a 64channel graph for each trial? - - % Plot the data - figure; - num_rows = size(event_triggered_lfps,2); - LFPs_per_shank = num_rows/ 8; % will be 8 for 64 channels, 7 for 56 channels (diff) - for i_row = 1 : num_rows - - y_lim = [-1500, 1500]; - plot_col = ceil(i_row / LFPs_per_shank); - plot_row = i_row - LFPs_per_shank * (plot_col-1); - plot_num = (plot_row-1) * 8 + plot_col; - - subplot(LFPs_per_shank,8,plot_num); - - % sue the trials structure data - % (trials_validchannels_marked.is_channel_valid_ordered) - % to color the actual plots based on good or - % bad sites above the specified threshold (1500mV, check the script to verify) -% switch trials_validchannels_marked(i_trial).is_channel_valid_ordered(i_row) - switch trials_validchannels_marked(trial_idx).is_channel_valid_ordered(i_row) - % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps - case 0 - plot_color = 'r'; % marks bad channels within specified trial - case 1 - plot_color = 'k'; % marks good channels within specified trial - otherwise - plot_color = 'b'; % catch in case the data was not input into the structure - end - - plot_channel_lfps = plot(t, channel_lfps(i_row, :), plot_color); % change to log10 -- plot(f, 10*log10(power_lfps(:,1))) - set(gca, 'ylim',y_lim); - grid on - - % give titles to each subplot - if contains(ratID, NN8x8) % if the ratID is in the list, it'll assign it the correct probe_type for ordering the LFP data correctly - caption = sprintf('NN8x8 #%d', site_order(i_row)); - elseif contains(ratID, ASSY156) - caption = sprintf('ASSY156 #%d', site_order(i_row)); - elseif contains(ratID, ASSY236) - caption = sprintf('ASSY236 #%d', site_order(i_row)); - end - title(caption, 'FontSize', 8); - - if plot_row < LFPs_per_shank - set(gca,'xticklabels',[]) - end - - if plot_col > 1 - set(gca,'yticklabels',[]) - end - - ax = gca; - % This section is coded to color the axes of - % the plots when checking the amplifier.dat - % files 'by eye' using Neuroscope - switch valid_sites_reordered(i_row) % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps - case 0 - ax.XColor = 'r'; % Red % marks bad channels within specified trial - ax.YColor = 'r'; % Red - % ax.ylabel = 'k'; - case 1 - ax.XColor = 'k'; % black % marks good channels within specified trial - ax.YColor = 'k'; % black - % ax.ylabel = 'k'; - case 2 - ax.XColor = 'b'; % blue % marks channels as 'variable' and could be good for portions of the whole amplifier.dat file but bad for others. Thus some channels may be good for only some trials, not all. - ax.YColor = 'b'; - % ax.ylabel = 'k'; - otherwise - ax.XColor = 'b'; % blue % catch in case the data was not input into the structure - ax.YColor = 'b'; - % ax.ylabel = 'k'; - end - end - set(gcf, 'WindowState', 'maximize'); % maximizes the window so that it exports the graphics with appropriate font size - - % Gives the whole plot page a name - A=cell(1,5); - A{1} = ['Subject: ' ratID]; - A{2} = ['Session: ' session_name]; - A{3} = ['Task Level: ' choiceRTdifficulty{logData.taskLevel+1}]; - A{4} = ['Trial Number: ' num2str(trial_idx)]; - A{5} = ['EventFieldname and Trial Type: ' eventFieldnames{i_event} '_' trialType]; - - sgtitle(A, 'Interpreter','none'); - - if i_trial == 1 - exportgraphics(gcf, trials_plot_full_name); - else - exportgraphics(gcf, trials_plot_full_name, 'append', true); %ADD APPEND FLAG HERE - end - close; - end - end - end -end \ No newline at end of file + end +end diff --git a/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots_TESTING.m b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots_TESTING.m new file mode 100644 index 00000000..c774079b --- /dev/null +++ b/LFPs/Analysis_code_JM/script_generate_trials_structure_bad_data_plots_TESTING.m @@ -0,0 +1,405 @@ +% script to write field potentials around each trial + +% probe_mapping_fname = '/Volumes/SharedX/Neuro-Leventhal/data/ChoiceTask/Probe Histology Summary/ProbeSite_Mapping.xlsx'; + +intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +rats_with_intan_sessions = find_rawdata_folders(intan_parent_directory); +valid_trials_original_folder = find_trials_struct_original_folders(intan_parent_directory); + +%% +fname = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary\Rat_Information_channels_to_discard.xlsx'; % for channels to ignore based on visualizing in Neuroscope etc +num_trials_to_plot = 5; + +% start with an LFP file +% get the trial structure for that session +% run this analysis + +trialType = ('allgo'); % to pull out trIdx of the trials structure (all trials) + +switch lower(trialType) + case 'allgo' + eventFieldnames = {'centerIn'}; + case 'correctgo' + eventFieldnames = {'cueOn', 'centerIn', 'centerOut', 'tone', 'sideIn', 'sideOut', 'foodClick', 'foodRetrievel'}; +end + + % eventlist = {'centerin','cueon','centerout', 'sidein', 'sideout', + % 'foodretrievel'}; This is the full eventlist for a correctGo trial. + % Choose one for generating event_triggered_lfps. +num_events = length(eventFieldnames); +t_win = [-2.5, 5]; % need this line for the event_triggered_lfps to select the correct + + +% lists for ratID probe_type +NN8x8 = ["R0326", "R0327", "R0372", "R0379", "R0374", "R0376", "R0378", "R0394", "R0395", "R0396", "R0412", "R0413"]; % Specify list of ratID associated with each probe_type +ASSY156 = ["R0411", "R0419"]; +ASSY236 = ["R0420", "R0425", "R0427", "R0457"]; + +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a', 'R0427_20220920a'}; +sessions_to_ignore1 = {'R0425_20220728_ChVE_220728_112601', 'R0427_20220920_Testing_220920_150255'}; +sessions_to_ignore2 = {'R0427_20220908a', 'R0427_20220909a', 'R0427_20220912a','R0427_20220913a', 'R0427_20220914a', 'R0427_20220915a', 'R0427_20220916a'}; +% Trying this as a workaround. Code wouldn't skip these two trials. R0425 - 15 hour session and R0427 no data (files didn't save correctly)? + +naming_convention; % for labeling graphs + +% choiceTask difficulty levels +choiceRTdifficulty = cell(1, 10); +choiceRTdifficulty{1} = 'poke any'; +choiceRTdifficulty{2} = 'very easy'; +choiceRTdifficulty{3} = 'easy'; +choiceRTdifficulty{4} = 'standard'; +choiceRTdifficulty{5} = 'advanced'; +choiceRTdifficulty{6} = 'choice VE'; +choiceRTdifficulty{7} = 'choice easy'; +choiceRTdifficulty{8} = 'choice standard'; +choiceRTdifficulty{9} = 'choice advanced'; +choiceRTdifficulty{10} = 'testing'; + +%% + %Testing just loading in the trials structure here + +for i_ratfolder = 1 : length(rats_with_intan_sessions) + + session_trials_struct_folders = valid_trials_original_folder(i_ratfolder).trials_folders; + intan_folders = rats_with_intan_sessions(i_ratfolder).intan_folders; + + for i_sessionfolder = 1 : length(intan_folders) + % extract the ratID and session name from the LFP file + session_path = intan_folders{i_sessionfolder}; + rd_metadata = parse_rawdata_folder(intan_folders{i_sessionfolder}); + session_trials_folder = create_trials_structure_data_folder(rd_metadata, intan_parent_directory); + ratID = rd_metadata.ratID; + pd_trials_data = parse_trials_struct_folder(session_path); + intan_session_name = rd_metadata.session_name; + session_name = rd_metadata.session_name; + + if any(strcmp(session_path, sessions_to_ignore)) % can't quite get this to debug but seems ok - it keeps running these sessions and catching errors (hence the need to skip them!) + continue; + end + +% if contains(ratID, NN8x8) || contains(ratID, ASSY156) % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes +% continue; +% end + +% if contains(ratID, 'R0326') || contains(ratID, 'R0372') || contains(ratID, 'R0376')... +% || contains(ratID, 'R0374') || contains(ratID, 'R0378') || contains(ratID, 'R0379') % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes +% continue; +% end + + +% if contains(ratID, NN8x8)|| contains(ratID, ASSY156)|| contains(ratID, 'R0420')|| contains(ratID, 'R0425') +% continue; +% end + + if contains(ratID, 'R0328') || contains(ratID, 'R0327') || contains(ratID, 'R0411') || contains(ratID, 'R0456') % the first style it wouldn't skip these sessions so trying it as the 'intan' name instead of just the rawdata folder name. + continue; % R0328 has no actual ephys; using these lines to skip unneeded data. R0327 Can't create trials struct; R0420 I haven't added lines for + end + + if contains(session_name, sessions_to_ignore) || contains(session_name, sessions_to_ignore1) || contains(ratID, 'DigiInputTest') % Just always ignore these sessions. R0411 no data, DigitInputTest is t est files + continue; + end + + % Load trials structure + for i_trialsfolder = 1:length(session_trials_struct_folders) + session_trials_structure = fullfile(intan_parent_directory, ratID, [ratID '-LFP-trials-structures-original'], session_name); + session_trials = find_trials_mat(session_trials_structure); + trials_structure_file = dir(fullfile(session_trials)); + trials_structure_fname = fullfile(trials_structure_file.folder, trials_structure_file.name); + trials_structure = load(trials_structure_fname); + trials = trials_structure.trials; + end + parentFolder = fullfile(intan_parent_directory, ... + ratID, ... + [ratID '-processed']); + + if contains(ratID, NN8x8) % if the ratID is in the list, it'll assign it the correct probe_type for ordering the LFP data correctly + probe_type = 'NN8x8'; + elseif contains(ratID, ASSY156) + probe_type = 'ASSY156'; + elseif contains(ratID, ASSY236) + probe_type = 'ASSY236'; + end + + % load_channel_information - this file is coded 0 = bad, 1 = good, + % 2 = variable for data in each channel for each session_name for + % each rat_ID. Use the opts.VariableNamesRange for eat ratID to + % detectImportOptions otherwise there's an error due to different + % session number for each rat + sheetname = ratID; + if contains(ratID, 'R0326') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:H1', 'datarange', 'A2:H65', 'sheet', sheetname); + elseif contains(ratID, 'R0327') || contains(ratID, 'R0374') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:E1', 'datarange', 'A2:E65', 'sheet', sheetname); + elseif contains(ratID, 'R0372') || contains(ratID, 'R0378')|| contains(ratID, 'R0396') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:J1', 'datarange', 'A2:J65', 'sheet', sheetname); + elseif contains(ratID, 'R0379') || contains(ratID, 'R0413') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:L1', 'datarange', 'A2:L65', 'sheet', sheetname); + elseif contains(ratID, 'R0376') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:O1', 'datarange', 'A2:O65', 'sheet', sheetname); + elseif contains(ratID, 'R0394') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:G1', 'datarange', 'A2:G65', 'sheet', sheetname); + elseif contains(ratID, 'R0395') || contains(ratID, 'R0427') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:K1', 'datarange', 'A2:K65', 'sheet', sheetname); + elseif contains(ratID, 'R0412') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:M1', 'datarange', 'A2:M65', 'sheet', sheetname); + elseif contains(ratID, 'R0419') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:P1', 'datarange', 'A2:P65', 'sheet', sheetname); + elseif contains(ratID, 'R0420') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:N1', 'datarange', 'A2:N65', 'sheet', sheetname); + elseif contains(ratID, 'R0425') + opts = detectImportOptions(fname, 'filetype', 'spreadsheet', 'VariableNamesRange', 'A1:V1', 'datarange', 'A2:V65', 'sheet', sheetname); + end + + probe_channel_info = load_channel_information(fname, sheetname); + [channel_information, intan_site_order, site_order] = channel_by_probe_site_ALL(probe_channel_info, probe_type); + +% trials_structure = [parentFolder(1:end-9) 'LFP-trials-structures']; +% if ~exist(trials_structure, 'dir') +% mkdir(trials_structure); +% end + + % Graphs in this folder are WITH calculating an outlier threshold. + trials_structure_plots = [parentFolder(1:end-9) 'LFP-trials-structures-graphs']; + if ~exist(trials_structure_plots, 'dir') + mkdir(trials_structure_plots); + end + + % Graphs in this folder are grandfathered in BEFORE I calculated a specific outlier threshold and just using by eye on Neuroscope. + processed_graphFolder_LFP_eventFieldname = [parentFolder(1:end-9) 'LFP-eventFieldname-graphs']; + if ~exist(processed_graphFolder_LFP_eventFieldname, 'dir') + mkdir(processed_graphFolder_LFP_eventFieldname); + end + +% % Start of generation of trials structure + [session_folder, ~, ~] = fileparts(intan_folders{i_sessionfolder}); % even with loading in trials structure, still need logData but it is quick to load. + session_log = find_session_log(session_folder); + + if isempty(session_log) + sprintf('no log file found for %s', session_folder) + end + + logData = readLogData(session_log); %gathering logData information +% +% % calculate nexData, need digital input and analog input files +% digin_fname = fullfile(intan_folders{i_sessionfolder}, 'digitalin.dat'); +% analogin_fname = fullfile(intan_folders{i_sessionfolder}, 'analogin.dat'); +% rhd_fname = fullfile(intan_folders{i_sessionfolder}, 'info.rhd'); +% +% if ~exist(digin_fname, 'file') +% sprintf('no digital input file for %s', session_folder); +% continue +% end +% +% if ~exist(analogin_fname, 'file') +% sprintf('no analog input file for %s', session_folder); +% continue +% end +% +% if ~exist(rhd_fname, 'file') +% sprintf('no rhd info file for %s', session_folder); +% continue +% end +% +% % read in rhd info; requires 'info.rhd' file. +% rhd_info = read_Intan_RHD2000_file_DL(rhd_fname); +% +% % read digital input file +% dig_data = readIntanDigitalFile(digin_fname); +% if ~isfield(rhd_info, 'board_adc_channels') +% sprintf('board_adc_channels field not found in rhd_info for %s', session_folder) +% continue +% end +% +% analog_data = readIntanAnalogFile(analogin_fname, rhd_info.board_adc_channels); +% nexData = intan2nex(dig_data, analog_data, rhd_info); +% +% try +% sprintf('attempting to create trials structure for %s', session_folder) +% trials = createTrialsStruct_simpleChoice_Intan( logData, nexData ); +% catch +% sprintf('could not generate trials structure for %s', session_folder) +% continue +% end +% % End of calculating trials structure (original without threshold info) + + + + % extract Trial Type (e.g. correctGo) + % Set trialType at beginning of file + [ trialEventParams ] = getTrialEventParams(trialType); % Need to verify this section. Working Here + trIdx = extractTrials(trials, trialEventParams); + num_trials = length(trIdx); + + %Getting the LFP-fname + pd_folder = create_processed_data_folder(rd_metadata, intan_parent_directory); + + lfp_fname = fullfile(pd_folder, create_lfp_fname(rd_metadata)); + lfp_data = load(lfp_fname); % Load in the LFP data for rearranging LFP data (ordering it by probe_type) + Fs = lfp_data.actual_Fs; % need the Fs loaded in for gathering event_triggered_lfps in the next 'for' loop + + % for now, skipping R0378_20210507a because the session only recorded 63 channels instead of 64. + % sessions_to_ignore doesn't seem to be working here so using this + % catch instead + lfp_data_num_rows = size(lfp_data.lfp,1); + if lfp_data_num_rows < 64 + continue; + end + + % Order the lfps here + [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); + % Orders the lfps by probe site mapping (double check the single file to remove catches for loading in single data) + + %find the column header that matches session_name (session_name and + %the column headers should match for each session to pull in the + %info for probe site health + valid_sites_reordered = channel_information.(session_name); + %channel_descriptions = find(probe_channel_info.Properties.VariableNames, session_name); + + % Generate event triggered LFPs + for i_event = 1 : num_events + % write code here to get event_triggered_lfps_ordered + % for each separate event + + % % create a filename to save the trials_structure + trials_fname_to_save = char(strcat(session_name, '_', eventFieldnames{i_event}, '_', trialType, '_', 'trials', '.mat')); % add in ''_', sprintf('trial%u.pdf', trial_idx)' should you want to save files individually + trials_full_name = fullfile(session_trials_folder, trials_fname_to_save); + + + % Cath for troubleshooting the data to NOT run through all of the folders and + % create the data if the plots are already made. This will not create 5 pages of 5 different random trials. + % It will simply generate one trial per session. + + % this function adds 2 fields to the trials structure to identify bad sites + trial_ts = extract_trial_ts(trials(trIdx), eventFieldnames{i_event}); % This is actually pulled in from the extract_event_related_LFPs script so do we need it here? + + event_triggered_lfps = extract_event_related_LFPs(ordered_lfp, trials(trIdx), eventFieldnames{i_event}, 'fs', Fs, 'twin', t_win); % made this using ordered_lfp data + tpoints_per_event = size(event_triggered_lfps, 3); + t = linspace(t_win(1), t_win(2), tpoints_per_event); + + trials_validchannels_marked = identify_bad_data(lfp_data, trials(trIdx), intan_site_order_for_trials_struct, probe_type, trialEventParams, eventFieldnames); + save(trials_full_name, 'trials_validchannels_marked'); + + % % at this point, event_triggered_lfps_ordered is an m x n x p array where + % % m is the number of events (i.e., all the cueon OR nosein OR other... + % % events) in a single session; n is the number of channels in that session + % % (generally 64), and p is the number of time samples extracted around each + % % event (i.e., p = 2001 if we extract +/- 2 seconds around each event). + + % above here, the lfp file, trial structure never change + + % create a filename to save the plots + trials_plot_fname_to_save = char(strcat(session_name, '_', eventFieldnames{i_event}, '_', trialType, '_', 'trials', '.pdf')); % add in ''_', sprintf('trial%u.pdf', trial_idx)' should you want to save files individually + trials_plot_full_name = fullfile(trials_structure_plots, trials_plot_fname_to_save); + + % % pts_per_event = size(event_triggered_lfps,3); % I think these lines need to go into the for loop for whenever a set of data is loaded into the workspace. + % % num_channels = size(event_triggered_lfps,2); + % num_trials = size(event_triggered_lfps, 1); + % + % catch for num_trials_to_plot < num_trials (this is specifically for R0420_20220714 since it did only 3 trials so can't graph 5 random trials). + if num_trials < num_trials_to_plot + continue; + end + + % select trials for plotting at random + trial_idx_to_plot = randperm(num_trials, num_trials_to_plot); % create a catch if num_trials_to_plot < num_trials, continue. + + for i_trial = 1:num_trials_to_plot + trial_idx = trial_idx_to_plot(i_trial); + channel_lfps = squeeze(event_triggered_lfps(trial_idx,:, :)); + % trials, create a 64channel graph for each trial? + + % Plot the data + figure; + num_rows = size(event_triggered_lfps,2); + LFPs_per_shank = num_rows/ 8; % will be 8 for 64 channels, 7 for 56 channels (diff) + for i_row = 1 : num_rows + + y_lim = [-2500, 2500]; + plot_col = ceil(i_row / LFPs_per_shank); + plot_row = i_row - LFPs_per_shank * (plot_col-1); + plot_num = (plot_row-1) * 8 + plot_col; + + subplot(LFPs_per_shank,8,plot_num); + + % sue the trials structure data + % (trials_validchannels_marked.is_channel_valid_ordered) + % to color the actual plots based on good or + % bad sites above the specified threshold (1500mV, check the script to verify) + % switch trials_validchannels_marked(i_trial).is_channel_valid_ordered(i_row) + switch trials_validchannels_marked(trial_idx).is_channel_valid_ordered(i_row) + % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps + case 0 + plot_color = 'r'; % marks bad channels within specified trial + case 1 + plot_color = 'k'; % marks good channels within specified trial + otherwise + plot_color = 'b'; % catch in case the data was not input into the structure + end + + plot_channel_lfps = plot(t, channel_lfps(i_row, :), plot_color); % change to log10 -- plot(f, 10*log10(power_lfps(:,1))) + set(gca, 'ylim',y_lim); + grid on + + % give titles to each subplot + if contains(ratID, NN8x8) % if the ratID is in the list, it'll assign it the correct probe_type for ordering the LFP data correctly + caption = sprintf('NN8x8 #%d', site_order(i_row)); + elseif contains(ratID, ASSY156) + caption = sprintf('ASSY156 #%d', site_order(i_row)); + elseif contains(ratID, ASSY236) + caption = sprintf('ASSY236 #%d', site_order(i_row)); + end + title(caption, 'FontSize', 8); + + if plot_row < LFPs_per_shank + set(gca,'xticklabels',[]) + end + + if plot_col > 1 + set(gca,'yticklabels',[]) + end + + ax = gca; + % This section is coded to color the axes of + % the plots when checking the amplifier.dat + % files 'by eye' using Neuroscope + switch valid_sites_reordered(i_row) % make sure is_valid_lfp is a boolean with true if it's a good channel; make sure this is in the same order as channel_lfps + case 0 + ax.XColor = 'r'; % Red % marks bad channels within specified trial + ax.YColor = 'r'; % Red + % ax.ylabel = 'k'; + case 1 + ax.XColor = 'k'; % black % marks good channels within specified trial + ax.YColor = 'k'; % black + % ax.ylabel = 'k'; + case 2 + ax.XColor = 'b'; % blue % marks channels as 'variable' and could be good for portions of the whole amplifier.dat file but bad for others. Thus some channels may be good for only some trials, not all. + ax.YColor = 'b'; + % ax.ylabel = 'k'; + otherwise + ax.XColor = 'b'; % blue % catch in case the data was not input into the structure + ax.YColor = 'b'; + % ax.ylabel = 'k'; + end + end + set(gcf, 'WindowState', 'maximize'); % maximizes the window so that it exports the graphics with appropriate font size + + % Gives the whole plot page a name + A=cell(1,5); + A{1} = ['Subject: ' ratID]; + A{2} = ['Session: ' session_name]; + A{3} = ['Task Level: ' choiceRTdifficulty{logData.taskLevel+1}]; + A{4} = ['Trial Number: ' num2str(trial_idx)]; + A{5} = ['EventFieldname and Trial Type: ' eventFieldnames{i_event} '_' trialType]; + + sgtitle(A, 'Interpreter','none'); + + if i_trial == 1 + exportgraphics(gcf, trials_plot_full_name); + else + exportgraphics(gcf, trials_plot_full_name, 'append', true); %ADD APPEND FLAG HERE + end + close; + end + end + end +end diff --git a/LFPs/Analysis_code_JM/utils/create_trials_structure_original_folder.m b/LFPs/Analysis_code_JM/utils/create_trials_structure_original_folder.m new file mode 100644 index 00000000..fe8b9aae --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/create_trials_structure_original_folder.m @@ -0,0 +1,10 @@ +function session_trials_folder_original = create_trials_structure_original_folder(rd_metadata, parent_directory) + +trials_folder= strcat(rd_metadata.ratID, '-LFP-trials-structures-original'); +session_trials_folder_original = fullfile(parent_directory, rd_metadata.ratID, trials_folder, rd_metadata.session_name); + +if ~isfolder(session_trials_folder_original) + + mkdir(session_trials_folder_original) + +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/utils/find_trials_struct_folders.m b/LFPs/Analysis_code_JM/utils/find_trials_struct_folders.m index f0bf3110..313228ad 100644 --- a/LFPs/Analysis_code_JM/utils/find_trials_struct_folders.m +++ b/LFPs/Analysis_code_JM/utils/find_trials_struct_folders.m @@ -1,12 +1,12 @@ -function valid_trials_folder = find_trials_struct_folders(intan_choicetask_parent) +function valid_trials_folder = find_trials_struct_folders(intan_parent_directory) -potential_rat_folders = dir(intan_choicetask_parent); +potential_rat_folders = dir(intan_parent_directory); valid_trials_folder = struct('name',[], 'trials_structure_folders', []); num_valid_rat_folders = 0; for i_folder = 1 : length(potential_rat_folders) - full_path = fullfile(intan_choicetask_parent, potential_rat_folders(i_folder).name); + full_path = fullfile(intan_parent_directory, potential_rat_folders(i_folder).name); if ~isfolder(full_path) continue end diff --git a/LFPs/Analysis_code_JM/utils/find_trials_struct_original_folders.m b/LFPs/Analysis_code_JM/utils/find_trials_struct_original_folders.m new file mode 100644 index 00000000..037e8c8c --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/find_trials_struct_original_folders.m @@ -0,0 +1,79 @@ +function valid_trials_original_folder = find_trials_struct_original_folders(intan_parent_directory) + +potential_rat_folders = dir(intan_parent_directory); +valid_trials_original_folder = struct('name',[], 'trials_structure_folders', []); + +num_valid_rat_folders = 0; +for i_folder = 1 : length(potential_rat_folders) + + full_path = fullfile(intan_parent_directory, potential_rat_folders(i_folder).name); + if ~isfolder(full_path) + continue + end + % look for folders of format RXXXX + if isvalidratfolder(potential_rat_folders(i_folder).name) + num_valid_rat_folders = num_valid_rat_folders + 1; + rat_folders{num_valid_rat_folders} = full_path; + end + +end + +num_rat_folders_with_trials_data = 0; +for i_ratfolder = 1 : num_valid_rat_folders + + [root_path, cur_ratID, ext] = fileparts(rat_folders{i_ratfolder}); + + cur_trials_folder = fullfile(rat_folders{i_ratfolder}, strcat(cur_ratID, '-LFP-trials-structures-original')); + + if ~isfolder(cur_trials_folder) + continue + end + + potential_session_folders = dir(cur_trials_folder); + + % rewrite this section to find session folders with processed data + num_valid_sessionfolders = 0; + found_trials_data = false; + num_trials_sessionfolders = 0; + for i_sessionfolder = 1 : length(potential_session_folders) + if isvalidchoicesessionfolder(potential_session_folders(i_sessionfolder).name) + num_valid_sessionfolders = num_valid_sessionfolders + 1; + + % test if session folder contains lfp data + full_pd_path = fullfile(cur_trials_folder, potential_session_folders(i_sessionfolder).name); +% test_folders = dir(full_pd_path); +% for i_tf = 1 : length(test_folders) +% fp = fullfile(full_pd_path, test_folders(i_tf).name); +% +% if ~isfolder(fp) || length(test_folders(i_tf).name) < 5 +% continue +% end +% +% if ~isvalidratfolder(test_folders(i_tf).name(1:5)) +% continue +% end +% +% % if isbehavior_vi_folder(test_folders(i_tf).name) +% % continue +% % end + + trials_datafolder = is_intan_trials_structure_datafolder(full_pd_path); + if ~isempty(trials_datafolder) + num_trials_sessionfolders = num_trials_sessionfolders + 1; + trials_datafolders{num_trials_sessionfolders} = trials_datafolder; + found_trials_data = true; + end + +% end + end + end + + if found_trials_data + num_rat_folders_with_trials_data = num_rat_folders_with_trials_data + 1; + valid_trials_original_folder(num_rat_folders_with_trials_data).name = rat_folders{i_ratfolder}; + valid_trials_original_folder(num_rat_folders_with_trials_data).trials_folders = trials_datafolders; + end + +end + +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/utils/generate_trials_structure.m b/LFPs/Analysis_code_JM/utils/generate_trials_structure.m new file mode 100644 index 00000000..707b8c3e --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/generate_trials_structure.m @@ -0,0 +1,138 @@ +% script to generate trials structures +% This file needs to be tested 12/18/2022 + +intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +rats_with_intan_sessions = find_rawdata_folders(intan_parent_directory); + +%% +% lists for ratID probe_type +NN8x8 = ["R0326", "R0327", "R0372", "R0379", "R0374", "R0376", "R0378", "R0394", "R0395", "R0396", "R0412", "R0413"]; % Specify list of ratID associated with each probe_type +ASSY156 = ["R0411", "R0419"]; +ASSY236 = ["R0420", "R0425", "R0427", "R0457"]; + +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a', 'R0427_20220920a', 'R0427_20220919a'}; +sessions_to_ignore1 = {'R0425_20220728_ChVE_220728_112601', 'R0427_20220920_Testing_220920_150255'}; +sessions_to_ignore2 = {'R0427_20220908a', 'R0427_20220909a', 'R0427_20220912a','R0427_20220913a', 'R0427_20220914a', 'R0427_20220915a', 'R0427_20220916a'}; +% Trying this as a workaround. Code wouldn't skip these two trials. R0425 - 15 hour session and R0427 no data (files didn't save correctly)? + +%% +for i_rat = 1 : length(rats_with_intan_sessions) + + intan_folders = rats_with_intan_sessions(i_rat).intan_folders; + + for i_sessionfolder = 1 : length(intan_folders) + % extract the ratID and session name from the LFP file + session_path = intan_folders{i_sessionfolder}; + rd_metadata = parse_rawdata_folder(intan_folders{i_sessionfolder}); + session_trials_folder_original = create_trials_structure_original_folder(rd_metadata, intan_parent_directory); + ratID = rd_metadata.ratID; + intan_session_name = rd_metadata.session_name; + session_name = rd_metadata.session_name; + + if any(strcmp(session_path, sessions_to_ignore)) % can't quite get this to debug but seems ok - it keeps running these sessions and catching errors (hence the need to skip them!) + continue; + end + + +% if contains(ratID, 'R0326') || contains(ratID, 'R0372') || contains(ratID, 'R0376')... +% || contains(ratID, 'R0374') || contains(ratID, 'R0378') || contains(ratID, 'R0379') % just trying to skip some lines of data to get to the last set to debug. Uncomment out to run more trialTypes +% continue; +% end + + +% if contains(ratID, NN8x8)|| contains(ratID, ASSY156)|| contains(ratID, 'R0420')|| contains(ratID, 'R0425') +% continue; +% end + + if contains(ratID, 'R0328') || contains(ratID, 'R0327') || contains(ratID, 'R0411') % the first style it wouldn't skip these sessions so trying it as the 'intan' name instead of just the rawdata folder name. + continue; % R0328 has no actual ephys; using these lines to skip unneeded data. R0327 Can't create trials struct; R0420 I haven't added lines for + end + + if contains(session_name, sessions_to_ignore) || contains(intan_session_name, sessions_to_ignore1)|| contains(ratID, 'DigiInputTest') % Just always ignore these sessions. R0411 no data, DigitInputTest is t est files + continue; + end + + parentFolder = fullfile(intan_parent_directory, ... + ratID, ... + [ratID '-processed']); + +% if contains(ratID, NN8x8) % if the ratID is in the list, it'll assign it the correct probe_type for ordering the LFP data correctly +% probe_type = 'NN8x8'; +% elseif contains(ratID, ASSY156) +% probe_type = 'ASSY156'; +% elseif contains(ratID, ASSY236) +% probe_type = 'ASSY236'; +% end + + +% trials_structure = [parentFolder(1:end-9) 'LFP-trials-structures']; +% if ~exist(trials_structure, 'dir') +% mkdir(trials_structure); +% end + + trials_structure_original = [parentFolder(1:end-9) 'LFP-trials-structures-original']; + if ~exist(trials_structure_original, 'dir') + mkdir(trials_structure_original); + end + + [session_folder, ~, ~] = fileparts(intan_folders{i_sessionfolder}); + session_log = find_session_log(session_folder); + + trials_fname_to_save = char(strcat(session_name, '_', 'trials', '.mat')); % add in ''_', sprintf('trial%u.pdf', trial_idx)' should you want to save files individually + trials_original_full_name = fullfile(session_trials_folder_original, trials_fname_to_save); + + if exist(trials_original_full_name, 'file') + continue; + end + + if isempty(session_log) + sprintf('no log file found for %s', session_folder) + end + + logData = readLogData(session_log); %gathersing logData information + + % calculate nexData, need digital input and analog input files + digin_fname = fullfile(intan_folders{i_sessionfolder}, 'digitalin.dat'); + analogin_fname = fullfile(intan_folders{i_sessionfolder}, 'analogin.dat'); + rhd_fname = fullfile(intan_folders{i_sessionfolder}, 'info.rhd'); + + if ~exist(digin_fname, 'file') + sprintf('no digital input file for %s', session_folder); + continue + end + + if ~exist(analogin_fname, 'file') + sprintf('no analog input file for %s', session_folder); + continue + end + + if ~exist(rhd_fname, 'file') + sprintf('no rhd info file for %s', session_folder); + continue + end + + % read in rhd info; requires 'info.rhd' file. + rhd_info = read_Intan_RHD2000_file_DL(rhd_fname); + + % read digital input file + dig_data = readIntanDigitalFile(digin_fname); + if ~isfield(rhd_info, 'board_adc_channels') + sprintf('board_adc_channels field not found in rhd_info for %s', session_folder) + continue + end + + analog_data = readIntanAnalogFile(analogin_fname, rhd_info.board_adc_channels); + nexData = intan2nex(dig_data, analog_data, rhd_info); + + try + sprintf('attempting to create trials structure for %s', session_folder) + trials = createTrialsStruct_simpleChoice_Intan( logData, nexData ); + catch + sprintf('could not generate trials structure for %s', session_folder) + continue + end + + + save(trials_original_full_name, 'trials'); + end +end diff --git a/LFPs/Analysis_code_JM/utils/gude-2022-12-20.log b/LFPs/Analysis_code_JM/utils/gude-2022-12-20.log new file mode 100644 index 00000000..39f1c1ca --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/gude-2022-12-20.log @@ -0,0 +1,6 @@ +09:22:23:119 [ALWAYS] GUDE Logging Started +09:22:23:124 [INFO] SqliteResumeCache::SqliteResumeCache basedirectory is null/empty +09:22:23:125 [ERROR] SqliteResumeCache::CreateSqliteResumeCache creating resume cache pointer failure +09:22:23:126 [ERROR] Initialize resume cache pointer failure +09:22:23:127 [INFO] gude policy POLICY_USER limit (4,4,1) chunk 2097152, chunking (UL 0 DL 1), resXfer 1, adapt 0, NSURL 1, timeout 86400 +09:22:23:130 [INFO] Initialized -- gude-lib version: v0.12.1 app: gude diff --git a/LFPs/Analysis_code_JM/utils/identify_bad_data.m b/LFPs/Analysis_code_JM/utils/identify_bad_data.m index 6fd48a61..a5bc88f1 100644 --- a/LFPs/Analysis_code_JM/utils/identify_bad_data.m +++ b/LFPs/Analysis_code_JM/utils/identify_bad_data.m @@ -27,7 +27,7 @@ lfp = lfp_data.lfp; -outlier_thresh = 1500; % in mV +outlier_thresh = 2000; % in mV t_win = [-2.5 5]; % eventFieldnames = {'cueOn'}; diff --git a/LFPs/Analysis_code_JM/utils/load_trials_structure.m b/LFPs/Analysis_code_JM/utils/load_trials_structure.m index 97e6d1d9..05016875 100644 --- a/LFPs/Analysis_code_JM/utils/load_trials_structure.m +++ b/LFPs/Analysis_code_JM/utils/load_trials_structure.m @@ -1,31 +1,28 @@ % script to load in trials structure -intan_choicetask_parent = 'X:\Neuro-Leventhal\data\ChoiceTask'; +intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; -valid_trials_folder = find_trials_struct_folders(intan_choicetask_parent); % find the folders with the trials structures to lab each lfp with data calculated as 'bad'; % need to match this with the loaded in lfp data - -fname = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary\Rat_Information_channels_to_discard.xlsx'; % for channels to ignore based on visualizing in Neuroscope etc +valid_trials_original_folder = find_trials_struct_original_folders(intan_parent_directory); % find the folders with the trials structures to lab each lfp with data calculated as 'bad'; % need to match this with the loaded in lfp data %% -for i_ratfolder = 1 : length(valid_trials_folder) + for i_ratfolder = 1 : length(valid_trials_original_folder) + + session_trials_struct_folders = valid_trials_original_folder(i_ratfolder).trials_folders; + + for i_trials_folder = 1 : length(session_trials_struct_folders) + % extract the ratID and session name from the LFP file + session_path = session_trials_struct_folders{i_trials_folder}; + pd_trials_data = parse_trials_struct_folder(session_path); + ratID = pd_trials_data.ratID; + session_name = pd_trials_data.session_name; - session_trials_struct_folders = valid_trials_folder(i_ratfolder).trials_folders; - - for i_sessionfolder = 1 : length(session_trials_struct_folders) - % extract the ratID and session name from the LFP file - session_path = session_trials_struct_folders{i_sessionfolder}; - pd_trials_data = parse_trials_struct_folder(session_path); - ratID = pd_trials_data.ratID; - session_name = pd_trials_data.session_name; - - % Load trials structure - for i_trialsfolder = 1:length(session_trials_struct_folders) - session_trials_structure = fullfile(intan_choicetask_parent, ratID, [ratID '-LFP-trials-structures'], session_name); - session_trials = find_trials_mat(session_trials_structure); - trials_structure_file = dir(fullfile(session_trials)); - trials_structure_fname = fullfile(trials_structure_file.folder, trials_structure_file.name); - trials_structure = load(trials_structure_fname); - trials_validchannels_marked = trials_structure.trials_validchannels_marked; + % Load trials structure + for i_trialsfolder = 1:length(session_trials_struct_folders) + session_trials_structure = fullfile(intan_parent_directory, ratID, [ratID '-LFP-trials-structures-original'], session_name); + session_trials = find_trials_mat(session_trials_structure); + trials_structure_file = dir(fullfile(session_trials)); + trials_structure_fname = fullfile(trials_structure_file.folder, trials_structure_file.name); + trials_structure = load(trials_structure_fname); + trials = trials_structure.trials; + end + end end - - end -end diff --git a/LFPs/Analysis_code_JM/utils/load_trials_structure_identify_bad_data.m b/LFPs/Analysis_code_JM/utils/load_trials_structure_identify_bad_data.m new file mode 100644 index 00000000..059ecb7b --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/load_trials_structure_identify_bad_data.m @@ -0,0 +1,30 @@ +% script to load in trials structure +intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; + +valid_trials_folder = find_trials_struct_folders(intan_parent_directory); % find the folders with the trials structures to lab each lfp with data calculated as 'bad'; % need to match this with the loaded in lfp data + +fname = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary\Rat_Information_channels_to_discard.xlsx'; % for channels to ignore based on visualizing in Neuroscope etc + +%% + for i_ratfolder = 1 : length(valid_trials_folder) + + session_trials_struct_folders = valid_trials_folder(i_ratfolder).trials_folders; + + for i_trials_folder = 1 : length(session_trials_struct_folders) + % extract the ratID and session name from the LFP file + session_path = session_trials_struct_folders{i_trials_folder}; + pd_trials_data = parse_trials_struct_folder(session_path); + ratID = pd_trials_data.ratID; + session_name = pd_trials_data.session_name; + + % Load trials structure + for i_trialsfolder = 1:length(session_trials_struct_folders) + session_trials_structure = fullfile(intan_parent_directory, ratID, [ratID '-LFP-trials-structures'], session_name); + session_trials = find_trials_mat(session_trials_structure); + trials_structure_file = dir(fullfile(session_trials)); + trials_structure_fname = fullfile(trials_structure_file.folder, trials_structure_file.name); + trials_structure = load(trials_structure_fname); + trials_validchannels_marked = trials_structure.trials_validchannels_marked; + end + end + end diff --git a/LFPs/Analysis_code_JM/utils/probe_site_mapping_all_probes_DL.m b/LFPs/Analysis_code_JM/utils/probe_site_mapping_all_probes_DL.m new file mode 100644 index 00000000..8a520562 --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/probe_site_mapping_all_probes_DL.m @@ -0,0 +1,159 @@ +function intan_to_site_map_DL = probe_site_mapping_all_probes_DL(probe_type) + +% INPUTS: + % probe_type - this current mapping is for the NeuroNexus H64LP A8x8 probe. +% OUTPUTS: + % intan_to_site_map - structure containing reorganization of the lfp_amplifier + % variable ordering each shank from ventral to dorsal + +if strcmpi(probe_type, 'nn8x8') + %this is how the lfp.mat file is organized + % lfp_organization = [1:64]; + + %shank number - organized based on lfp.mat files generated in -processed + %folders to match the intan_ampflier number; remember Intan_amplifer + %numbers are 0:63, while matlab (aka .mat files) are 1:64 + % DL - this maps the numbering on the NeuroNexus omnetics connector to + % the numbering on the Intan amplifier + + % NOTE: site 35 is entered twice in the matrix below + intan_amplifier = [49:56,...%shank1, nn sites 1-8 + 57:64,...%shank2 + 33:35, 37:39, 41:42,...%shank3 + 35,40,43:48,...%shank4 + 17:22,26,30,...%shank5 + 23:25, 27:29, 31:32,...%shank6 + 1:8,...%shank7 + 9:16]; %shank8 + + intan_amplifier_DL = [49, 48, 51, 50, 53, 52, 55, 54, ... + 57, 56, 59, 58, 61, 60, 63, 62, ... + 32, 33, 34, 36, 37, 38, 40, 41, ... + 42, 44, 45, 46, 47, 43, 39, 35, ... + 29, 25, 21, 17, 16, 19, 18, 20, ... + 23, 22, 24, 27, 26, 28, 31, 30, ... + 00, 01, 02, 03, 04, 05, 06, 07, ... + 08, 09, 10, 11, 12, 13, 14, 15]; + % above is the actual intan channel numbers. Add one because Matlab + % starts indexing at 1. Note that if this were ported to Python, would + % NOT want to add 1 here + intan_amplifier_DL = intan_amplifier_DL + 1; + + +% Map sites ventral to dorsal on the NeuroNexus H64LP A8x8 probe. +% THIS SECTION ORDERS THE NN-SITE BASED ON THE ACTUAL AMPLIFIER ORDER FROM +% THE LFP.MAT FILE OR REARRANGE FROM LINE 14-21 ABOVE. +% DL - this maps the ordering on the NN shank to the NN omnetics connector. +% Therefore, intan_amplifier(NNsite_order) gives the row in the original +% LFP file that corresponds to spatially ordered sites from VENTRAL to +% DORSAL + NNsite_order = [2,7,1,8,4,5,3,6,... %shank 1 + 10,15,9,16,12,13,11,14,... %shank2 + 17,24,18,23,19,22,20,21,...%shank3 + 27,25,29,26,30,28,31,32,...%shank4 + 40,37,39,35,38,36,34,33,...%shank5 + 42,47,41,48,43,46,45,44,...%shank6 + 49,56,50,55,51,54,52,53,...%shank7 + 57,64,58,63,59,62,60,61]; %shank8 + + NNsite_order_DL = [05, 04, 06, 03, 07, 02, 08, 01, ... + 13, 12, 14, 11, 15, 10, 16, 09, ... + 21, 20, 22, 19, 23, 18, 24, 17, ... + 29, 28, 30, 27, 31, 26, 32, 25, ... + 37, 36, 38, 35, 39, 34, 40, 33, ... + 45, 44, 46, 43, 47, 42, 48, 41, ... + 53, 52, 54, 51, 55, 50, 56, 49, ... + 61, 60, 62, 59, 63, 58, 64, 57]; + +% % Map sites ventral to dorsal on the NeuroNexus H64LP A8x8 probe. +% THIS CHUNK OF CODE IS THE LITERAL NEURONEXUS NUMBERS VENTRAL TO DORSAL ORDER. +% NNsite_order = [1,8,2,7,3,6,4,5,... %shank 1 +% 9,16,10,15,11,14,12,13,... %shank2 +% 17,24,18,23,19,22,20,21,...%shank3 +% 25,32,26,31,27,30,28,29,...%shank4 +% 33,40,34,39,35,38,36,37,...%shank5 +% 41,48,42,47,43,46,44,45,...%shank6 +% 49,56,50,55,51,54,52,53,...%shank7 +% 57,64,58,63,59,62,60,61];%shank8 + + intan_to_site_map = intan_amplifier(NNsite_order)'; % THIS IS VENTRAL TO DORSAL + intan_to_site_map_DL = intan_amplifier_DL(NNsite_order_DL)'; % THIS IS DORSAL TO VENTRAL + +elseif strcmpi(probe_type, 'assy156') + % Per Jen, Intan chip is facing down on the ASSY-156 probe, but this + % numbering suggests the chip was facing up. This doesn't seem to match + % with the way she did the NN mapping; for example, site 1 from + % Cambridge DOES NOT map to site 0 intan. I redid the mapping and they + % seem to match other than the duplication of site 36 in + % Cambridge156_order + intan_amplifier = [1:17,19:32,63,...% Shank A + 18,33:62, 64]; % Shank B + intan_amplifier_DL = [47, 46, 45, 44, 43, 42, 41, 40, ... + 39, 38, 37, 36, 35, 34, 33, 32, ... + 31, 30, 29, 28, 27, 26, 25, 24, ... + 23, 22, 21, 20, 19, 18, 17, 16, ... + 15, 14, 13, 12, 11, 10, 09, 08, ... + 07, 06, 05, 04, 03, 02, 01, 00, ... + 63, 62, 61, 60, 59, 58, 57, 56, ... + 55, 54, 53, 52, 51, 50, 49, 48]; + intan_amplifier_DL = intan_amplifier_DL + 1; + + Cambridge156_order_DL = [26, 35, 33, 24, 37, 28, 41, 39, ... + 27, 40, 29, 50, 47, 32, 19, 44, ... + 23, 36, 38, 25, 34, 21, 22, 43, ... + 30, 20, 45, 17, 48, 18, 46, 42, ... + 12, 53, 04, 14, 51, 03, 15, 01, ... + 31, 49, 06, 57, 58, 62, 07, 11, ... + 16, 02, 52, 13, 56, 54, 08, 61, ... + 63, 10, 59, 55, 05, 09, 60, 64]; + % note site 36 is entered twice in the matrix below +% Cambridge156_order = [22,14,16,24,12,20,8,10,...% Shank A. This is verified JM 11/5/2022 +% 21,9,19,32,2,17,29,5,... +% 25,13,11,23,36,27,26,6,... +% 18,28,4,31,1,30,3,7,... +% 38,61,46,36,63,47,35,49,...% Shank B +% 33,64,44,57,56,52,43,39,... +% 34,48,62,37,58,60,42,53,... +% 51,40,55,59,45,41,54,50]; + + % THIS GOES FROM DORSAL TO VENTRAL - DL +% intan_to_site_map = intan_amplifier(Cambridge156_order)'; + intan_to_site_map_DL = intan_amplifier_DL(Cambridge156_order_DL)'; + +elseif strcmpi(probe_type, 'assy236') + + % see Cambridge Neurotech Mini-Amp-64 User Guide + intan_amplifier_DL = [41, 39, 38, 37, 35, 34, 33, 32, ... + 29, 30, 28, 26, 25, 24, 22, 20, ... + 46, 45, 44, 43, 42, 40, 36, 31, ... + 27, 23, 21, 18, 19, 17, 16, 14, ... + 55, 53, 54, 52, 51, 50, 49, 48, ... + 47, 15, 13, 12, 11, 09, 10, 08, ... + 63, 62, 61, 60, 59, 58, 57, 56, ... + 07, 06, 05, 04, 03, 02, 01, 00]; + intan_amplifier_DL= intan_amplifier_DL + 1; + + % verified order below + Cambridge236_order_DL = [64, 48, 16, 61, 13, 60, 29, 12, ... + 59, 32, 42, 24, 09, 57, 44, 45, ... + 47, 14, 31, 62, 15, 63, 46, 30, ... + 58, 27, 10, 26, 25, 43, 28, 11, ... + 34, 18, 54, 23, 07, 38, 22, 40, ... + 56, 08, 52, 05, 17, 03, 51, 50, ... + 41, 55, 20, 36, 06, 39, 49, 19, ... + 01, 35, 04, 21, 53, 33, 37, 02]; +% intan_amplifier = [1:32,...% Shank A % Verified 11/5/2022 JM +% 33:64]; % Shank B +% Cambridge236_order = [1,9,21,4,26,5,20,27,...% Shank A % Verified 11/5/2022 +% 6,15,16,32,30,8,13,12,... +% 11,25,17,3,23,2,10,18,... +% 7,22,31,24,28,14,19,29,... +% 54,46,59,37,34,51,41,49,...%shankB +% 57,33,61,36,47,39,62,63,... +% 48,58,44,53,35,50,64,45,... +% 42,55,38,43,60,56,52,40]; + +% intan_to_site_map = intan_amplifier(Cambridge236_order)'; + intan_to_site_map_DL = intan_amplifier_DL(Cambridge236_order_DL)'; +end +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/utils/ratID_by_probe_type.xlsx b/LFPs/Analysis_code_JM/utils/ratID_by_probe_type.xlsx index 3bcc50b7..560c51d4 100644 Binary files a/LFPs/Analysis_code_JM/utils/ratID_by_probe_type.xlsx and b/LFPs/Analysis_code_JM/utils/ratID_by_probe_type.xlsx differ diff --git a/LFPs/Analysis_code_JM/utils/script_calculate_perievent_scalograms.m b/LFPs/Analysis_code_JM/utils/script_calculate_perievent_scalograms.m new file mode 100644 index 00000000..57b4cc12 --- /dev/null +++ b/LFPs/Analysis_code_JM/utils/script_calculate_perievent_scalograms.m @@ -0,0 +1,69 @@ +% script_calculate_perievent_scalograms + +% probe_mapping_fname = '/Volumes/SharedX/Neuro-Leventhal/data/ChoiceTask/Probe Histology Summary/ProbeSite_Mapping.xlsx'; + +intan_parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; + +rats_with_intan_sessions = find_rawdata_folders(intan_parent_directory); + +%% +for i_rat = 1 : length(rats_with_intan_sessions) + + intan_folders = rats_with_intan_sessions(i_rat).intan_folders; + + for i_sessionfolder = 1 : length(intan_folders) + + [session_folder, ~, ~] = fileparts(intan_folders{i_sessionfolder}); + session_log = find_session_log(session_folder); + + if isempty(session_log) + sprintf('no log file found for %s', session_folder) + end + logData = readLogData(session_log); + + % calculate nexData, need digital input and analog input files + digin_fname = fullfile(intan_folders{i_sessionfolder}, 'digitalin.dat'); + analogin_fname = fullfile(intan_folders{i_sessionfolder}, 'analogin.dat'); + rhd_fname = fullfile(intan_folders{i_sessionfolder}, 'info.rhd'); + + if ~exist(digin_fname, 'file') + sprintf('no digital input file for %s', session_folder) + continue + end + + if ~exist(analogin_fname, 'file') + sprintf('no analog input file for %s', session_folder) + continue + end + + if ~exist(rhd_fname, 'file') + sprintf('no rhd info file for %s', session_folder) + continue + end + + % read in rhd info + rhd_info = read_Intan_RHD2000_file_DL(rhd_fname); + + % read digital input file + dig_data = readIntanDigitalFile(digin_fname); + if ~isfield(rhd_info, 'board_adc_channels') + sprintf('board_adc_channels field not found in rhd_info for %s', session_folder) + continue + end + + analog_data = readIntanAnalogFile(analogin_fname, rhd_info.board_adc_channels); + nexData = intan2nex(dig_data, analog_data, rhd_info); + + try + sprintf('attempting to create trials structure for %s', session_folder) + trials = createTrialsStruct_simpleChoice_Intan( logData, nexData ); + catch + sprintf('could not generate trials structure for %s', session_folder) + continue + end + % update here for next steps in calculating perievent scalograms + + sprintf('placeholder') + end + +end \ No newline at end of file diff --git a/LFPs/Analysis_code_JM/~$P_Analysis_Workflow.docx b/LFPs/Analysis_code_JM/~$P_Analysis_Workflow.docx deleted file mode 100644 index fbbe2f8f..00000000 Binary files a/LFPs/Analysis_code_JM/~$P_Analysis_Workflow.docx and /dev/null differ diff --git a/LFPs/calculate_bipolar_LFPs.m b/LFPs/calculate_bipolar_LFPs.m new file mode 100644 index 00000000..3b942c75 --- /dev/null +++ b/LFPs/calculate_bipolar_LFPs.m @@ -0,0 +1,40 @@ +function [bipolar_LFPs, probe_site_mapping] = calculate_bipolar_LFPs(lfp, probe_type) +%UNTITLED Summary of this function goes here +% Detailed explanation goes here + +% re-arrange the rows to reflect the order of LFPs +probe_site_mapping = probe_site_mapping_all_probes(probe_type); + +% first, rearrange field potentials so they are sorted from dorsal to +% ventral + +sorted_lfps = lfp(probe_site_mapping, :); +num_lfps = size(lfp, 1); % assume each row is a recording channel +num_pts = size(lfp, 2); + +switch lower(probe_type) + case 'nn8x8' + sites_per_column = 8; + case 'assy156' + sites_per_column = 16; + case 'assy236' + sites_per_column = 16; +end +num_columns = num_lfps / sites_per_column; +num_bipolar_lfps = num_lfps - num_columns; + +bipolar_LFPs = zeros(num_bipolar_lfps, num_pts); + +for i_sitecol = 1 : num_columns + + start_LFP_row = (i_sitecol - 1) * sites_per_column + 1; + end_LFP_row = i_sitecol * sites_per_column; + + start_bipolar_row = (i_sitecol-1) * (sites_per_column-1) + 1; + end_bipolar_row = i_sitecol * (sites_per_column-1); + + bipolar_LFPs(start_bipolar_row:end_bipolar_row, :) = diff(sorted_lfps(start_LFP_row:end_LFP_row, :), 1, 1); + +end + +end \ No newline at end of file diff --git a/LFPs/calculate_trial_spectrograms.m b/LFPs/calculate_trial_spectrograms.m new file mode 100644 index 00000000..77c437d4 --- /dev/null +++ b/LFPs/calculate_trial_spectrograms.m @@ -0,0 +1,134 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_types = {'monopolar', 'bipolar'}; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +t_window = [-2.5, 2.5]; +event_list = {'cueOn', 'centerIn', 'tone', 'centerOut' 'sideIn', 'sideOut', 'foodClick', 'foodRetrieval'}; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + % load trials structure + trials_name = sprintf('%s_trials.mat', session_name); + trials_name = fullfile(cur_dir, trials_name); + + if ~exist(trials_name) + sprintf('no trials structure found for %s', session_name) + continue + end + + load(trials_name) + + selected_trials = extract_trials_by_features(trials, trials_to_analyze); + if isempty(selected_trials) + sprintf('no %s trials found for %s', trials_to_analyze, session_name) + continue + end + + lfp_fname = strcat(session_name, '_lfp.mat'); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_fname) + continue + end + + lfp_data = load(lfp_fname); + + Fs = lfp_data.actual_Fs; + samp_window = round(t_window * Fs); + samples_per_event = range(samp_window) + 1; + fb = cwtfilterbank('signallength', samples_per_event, ... + 'samplingfrequency', Fs, ... + 'wavelet','amor',... + 'frequencylimits', [1, 100]); + + [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); + + for i_lfptype = 1 : length(lfp_types) + + lfp_type = lfp_types{i_lfptype}; + if strcmpi(lfp_type, 'bipolar') + ordered_lfp = diff_lfp_from_monopolar(ordered_lfp, probe_type); + end + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + num_channels = size(ordered_lfp, 1); + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s, %s', session_name, event_name, lfp_type) + + perievent_data = extract_perievent_data(ordered_lfp, selected_trials, event_list{4}, t_window, Fs); + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + % scalo_folder = sprintf('%s_scalos_%s', session_name, event_name); + % scalo_folder = fullfile(cur_dir, scalo_folder); + if ~exist(scalo_folder, 'dir') + mkdir(scalo_folder) + end + + for i_channel = 1 : num_channels + [shank_num, site_num] = get_shank_and_site_num(probe_lfp_type, i_channel); + + scalo_name = sprintf('%s_scalos_%s_%s_shank%02d_site%02d.mat',session_name, lfp_type, event_name, shank_num, site_num); + scalo_name = fullfile(scalo_folder, scalo_name); + + if exist(scalo_name, 'file') + continue + end + + event_triggered_lfps = squeeze(perievent_data(:, i_channel, :)); + + % comment back in if running on a machine without a gpu +% disp('cpu') +% tic +% [event_related_scalos, ~, coi] = trial_scalograms(event_triggered_lfps, fb); +% toc + + etl_g = gpuArray(event_triggered_lfps); + [event_related_scalos, ~, coi] = trial_scalograms(etl_g, fb); + + save(scalo_name, 'event_related_scalos', 'event_triggered_lfps', 'fb', 'coi', 't_window', 'i_channel'); + % saving i_channel is a check to make sure that shank + % and site are numbered correctly later + + end + end + end + end + +end \ No newline at end of file diff --git a/LFPs/calculate_trial_spectrograms_from_saved_bipolars.asv b/LFPs/calculate_trial_spectrograms_from_saved_bipolars.asv new file mode 100644 index 00000000..9386f40d --- /dev/null +++ b/LFPs/calculate_trial_spectrograms_from_saved_bipolars.asv @@ -0,0 +1,137 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'Z:\data\ChoiceTask\'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'Z:\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_types = {'monopolar', 'bipolar'}; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +t_window = [-2.5, 2.5]; +event_list = {'cueOn', 'centerIn', 'tone', 'centerOut' 'sideIn', 'sideOut', 'foodClick', 'foodRetrieval'}; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.RatID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + % load trials structure + trials_name = sprintf('%s_trials.mat', session_name); + trials_name = fullfile(cur_dir, trials_name); + + if ~exist(trials_name) + sprintf('no trials structure found for %s', session_name) + continue + end + + load(trials_name) + + selected_trials = extract_trials_by_features(trials, trials_to_analyze); + if isempty(selected_trials) + sprintf('no %s trials found for %s', trials_to_analyze, session_name) + continue + end + + lfp_fname = strcat(session_name, '_lfp.mat'); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_fname) + continue + end + + lfp_data = load(lfp_fname); + + Fs = lfp_data.actual_Fs; + samp_window = round(t_window * Fs); + samples_per_event = range(samp_window) + 1; + fb = cwtfilterbank('signallength', samples_per_event, ... + 'samplingfrequency', Fs, ... + 'wavelet','amor',... + 'frequencylimits', [1, 100]); + + [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); + + for i_lfptype = 1 : length(lfp_types) + + lfp_type = lfp_types{i_lfptype}; + if strcmpi(lfp_type, 'bipolar') + lfp_fname = strcat(session_name, '_bipolar_lfp.mat'); + else + lfp_fname = strcat(session_name, '_lfp.mat'); + end + lfp_fname = fullfile(cur_dir, ) + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + num_channels = size(ordered_lfp, 1); + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s, %s', session_name, event_name, lfp_type) + + perievent_data = extract_perievent_data(ordered_lfp, selected_trials, event_list{4}, t_window, Fs); + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + % scalo_folder = sprintf('%s_scalos_%s', session_name, event_name); + % scalo_folder = fullfile(cur_dir, scalo_folder); + if ~exist(scalo_folder, 'dir') + mkdir(scalo_folder) + end + + for i_channel = 1 : num_channels + [shank_num, site_num] = get_shank_and_site_num(probe_lfp_type, i_channel); + + scalo_name = sprintf('%s_scalos_%s_%s_shank%02d_site%02d.mat',session_name, lfp_type, event_name, shank_num, site_num); + scalo_name = fullfile(scalo_folder, scalo_name); + +% if exist(scalo_name, 'file') +% continue +% end + + event_triggered_lfps = squeeze(perievent_data(:, i_channel, :)); + + % comment back in if running on a machine without a gpu +% disp('cpu') +% tic +% [event_related_scalos, ~, coi] = trial_scalograms(event_triggered_lfps, fb); +% toc + + etl_g = gpuArray(event_triggered_lfps); + [event_related_scalos, ~, coi] = trial_scalograms(etl_g, fb); + + save(scalo_name, 'event_related_scalos', 'event_triggered_lfps', 'fb', 'coi', 't_window', 'i_channel'); + % saving i_channel is a check to make sure that shank + % and site are numbered correctly later + + end + end + end + end + +end \ No newline at end of file diff --git a/LFPs/calculate_trial_spectrograms_from_saved_bipolars.m b/LFPs/calculate_trial_spectrograms_from_saved_bipolars.m new file mode 100644 index 00000000..b3dfc059 --- /dev/null +++ b/LFPs/calculate_trial_spectrograms_from_saved_bipolars.m @@ -0,0 +1,151 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_types = {'monopolar', 'bipolar'}; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +t_window = [-2.5, 2.5]; +event_list = {'cueOn', 'centerIn', 'tone', 'centerOut' 'sideIn', 'sideOut', 'foodClick', 'foodRetrieval'}; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + % load trials structure + trials_name = sprintf('%s_trials.mat', session_name); + trials_name = fullfile(cur_dir, trials_name); + + if ~exist(trials_name) + sprintf('no trials structure found for %s', session_name) + continue + end + + load(trials_name) + + selected_trials = extract_trials_by_features(trials, trials_to_analyze); + if isempty(selected_trials) + sprintf('no %s trials found for %s', trials_to_analyze, session_name) + continue + end + + lfp_fname = strcat(session_name, '_lfp.mat'); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_fname) + continue + end + + lfp_data = load(lfp_fname); + + Fs = lfp_data.actual_Fs; + samp_window = round(t_window * Fs); + samples_per_event = range(samp_window) + 1; + fb = cwtfilterbank('signallength', samples_per_event, ... + 'samplingfrequency', Fs, ... + 'wavelet','amor',... + 'frequencylimits', [1, 100]); + +% [ordered_lfp, intan_site_order, intan_site_order_for_trials_struct, site_order] = lfp_by_probe_site_ALL(lfp_data, probe_type); + + for i_lfptype = 1 : length(lfp_types) + + lfp_type = lfp_types{i_lfptype}; + if strcmpi(lfp_type, 'bipolar') + lfp_fname = strcat(session_name, '_bipolar_lfp.mat'); + else + lfp_fname = strcat(session_name, '_lfp.mat'); + end + lfp_fname = fullfile(cur_dir, lfp_fname); + if ~isfile(lfp_fname) + sprintf('%s not found, skipping', lfp_fname) + continue + end + + lfp_data = load(lfp_fname); + + if strcmpi(lfp_type, 'bipolar') + ordered_lfp = lfp_data.bipolar_lfp; + else + probe_site_mapping = probe_site_mapping_all_probes(probe_type); + ordered_lfp = lfp_data.lfp(probe_site_mapping, :); + end + + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + num_channels = size(ordered_lfp, 1); + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s, %s', session_name, event_name, lfp_type) + + perievent_data = extract_perievent_data(ordered_lfp, selected_trials, event_name, t_window, Fs); + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + % scalo_folder = sprintf('%s_scalos_%s', session_name, event_name); + % scalo_folder = fullfile(cur_dir, scalo_folder); + if ~exist(scalo_folder, 'dir') + mkdir(scalo_folder) + end + + for i_channel = 1 : num_channels + [shank_num, site_num] = get_shank_and_site_num(probe_lfp_type, i_channel); + + scalo_name = sprintf('%s_scalos_%s_%s_shank%02d_site%02d.mat',session_name, lfp_type, event_name, shank_num, site_num); + scalo_name = fullfile(scalo_folder, scalo_name); + + if exist(scalo_name, 'file') + continue + end + + event_triggered_lfps = squeeze(perievent_data(:, i_channel, :)); + + % comment back in if running on a machine without a gpu + disp('cpu') + tic + [event_related_scalos, ~, coi] = trial_scalograms(event_triggered_lfps, fb); + toc + + etl_g = gpuArray(event_triggered_lfps); + [event_related_scalos, ~, coi] = trial_scalograms(etl_g, fb); + + save(scalo_name, 'event_related_scalos', 'event_triggered_lfps', 'fb', 'coi', 't_window', 'i_channel'); + % saving i_channel is a check to make sure that shank + % and site are numbered correctly later + + end + end + end + end + +end \ No newline at end of file diff --git a/LFPs/create_scalogram_map.m b/LFPs/create_scalogram_map.m new file mode 100644 index 00000000..392d7169 --- /dev/null +++ b/LFPs/create_scalogram_map.m @@ -0,0 +1,63 @@ +function [outputArg1,outputArg2] = create_scalogram_map(session_name,event_name,lfp_type,parent_directory) +%UNTITLED Summary of this function goes here +% function to create a map of scalograms for a specific event across all +% channels + +scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); + +h_fig = figure('papertype','usletter',... + 'PaperOrientation','landscape', ... + 'units', 'inches',... + 'position',[0 0 11 8.5]); +fig_title = sprintf('%s, %s, %s', session_name, event_name, lfp_type); + +num_cols = 8; % may need to change for the Cambridge probes +switch lfp_type + case 'monopolar' + num_rows = 8; + fname_append = ''; + case 'bipolar' + num_rows = 7; + fname_append = 'bipolar'; +end +t_layout = tiledlayout(num_rows, 8, 'tilespacing','tight'); + +for ii = 1 : 64 + shank_num = ceil(ii/num_rows); + site_num = ii - (shank_num-1) * num_rows; + + scalo_name = sprintf('%s_scalos_%s_shank%02d_site%02d.mat',session_name, event_name, shank_num, site_num); + scalo_name = fullfile(scalo_folder, scalo_name); + + load(scalo_name) +% % f = centerFrequencies(fb); +% + mean_scalo = squeeze(mean(log(abs(event_related_scalos)), 1)); + + tile_num = (site_num-1) * num_cols + shank_num; + ax = nexttile(tile_num); + display_scalogram(mean_scalo, t_window, fb, 'ax', ax); +% plot(sin(site_num*(1:10))) + if shank_num > 1 + set(gca,'yticklabel','') + end + if site_num < num_rows + set(gca,'XTickLabel','') + end +end + +title(t_layout, fig_title, 'interpreter', 'none') +xlabel(t_layout, 'time (s)') +ylabel(t_layout, 'frequency (Hz)') + +savename = sprintf('%s_%s_%s_meanlogpower.pdf', session_name, event_name, lfp_type); +fig_savename = sprintf('%s_%s_%s_meanlogpower.fig', session_name, event_name, lfp_type); +savename = fullfile(scalo_folder, savename); +fig_savename = fullfile(scalo_folder, fig_savename); + +print(savename, '-dpdf') +savefig(h_fig, fig_savename, 'compact') + +close(h_fig) + +end \ No newline at end of file diff --git a/LFPs/diff_lfp_from_monopolar.m b/LFPs/diff_lfp_from_monopolar.m new file mode 100644 index 00000000..735aeff2 --- /dev/null +++ b/LFPs/diff_lfp_from_monopolar.m @@ -0,0 +1,39 @@ +function diff_lfps = diff_lfp_from_monopolar(monopolar_ordered_lfp,probe_type) +%UNTITLED8 Summary of this function goes here +% INPUTS +% monopolar_ordered_lfp: num_channels x num_points array +% probe_type: string containing probe type. possibilities are 'nn8x8', + +num_sites = size(monopolar_ordered_lfp, 1); +num_points = size(monopolar_ordered_lfp, 2); +switch lower(probe_type) + case 'nn8x8' + % assume lfp is ordered such that monopolar_ordered_lfp(1:8,:) is + % the signal from dorsal to ventral on shank 1, + % monopolar_ordered_lfp(9:16,:) is the signal from dorsal to + % ventral on shank 2, etc. + + sites_per_shank = 8; + diff_sites_per_shank = sites_per_shank - 1; + num_shanks = num_sites / sites_per_shank; + + num_diff_sites = diff_sites_per_shank * num_shanks; + diff_lfps = zeros(num_diff_sites, num_points); + + for i_shank = 1 : num_shanks + diff_chan_start_idx = (i_shank-1) * diff_sites_per_shank + 1; + diff_chan_end_idx = i_shank * diff_sites_per_shank; + + mono_chan_start_idx = (i_shank-1) * sites_per_shank + 1; + mono_chan_end_idx = i_shank * sites_per_shank; + + % take differences going down columns + diff_lfps(diff_chan_start_idx:diff_chan_end_idx, :) = ... + diff(monopolar_ordered_lfp(mono_chan_start_idx:mono_chan_end_idx,:),1,1); + + end + +end + + +end \ No newline at end of file diff --git a/LFPs/script_display_session_trial_spectrograms.m b/LFPs/script_display_session_trial_spectrograms.m new file mode 100644 index 00000000..43138dcc --- /dev/null +++ b/LFPs/script_display_session_trial_spectrograms.m @@ -0,0 +1,62 @@ +% script to calculate scalograms for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'Z:\data\ChoiceTask\'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'Z:\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +% change the line below to allow looping through multiple trial types, +% extract left vs right, etc. +trials_to_analyze = 'correct'; +lfp_type = 'monopolar'; +% trials_to_analyze = 'all'; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.RatID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + probe_lfp_type = sprintf('%s_%s', probe_type, lfp_type); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + for i_event = 1 : length(event_list) + event_name = event_list{i_event}; + sprintf('working on session %s, event %s', session_name, event_name) + + perievent_data = extract_perievent_data(ordered_lfp, selected_trials, event_list{4}, t_window, Fs); + + scalo_folder = create_scalo_folder(session_name, event_name, parent_directory); +% scalo_folder = sprintf('%s_scalos_%s', session_name, event_name); +% scalo_folder = fullfile(cur_dir, scalo_folder); + if ~exist(scalo_folder, 'dir') + continue + end + create_scalogram_map(session_name, event_name, lfp_type, parent_directory); + + end + end + +end \ No newline at end of file diff --git a/LFPs/script_extract_LFPs_DL.asv b/LFPs/script_extract_LFPs_DL.asv new file mode 100644 index 00000000..9c51ca5e --- /dev/null +++ b/LFPs/script_extract_LFPs_DL.asv @@ -0,0 +1,72 @@ +% script to calculate LFPs for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a','R0427_20220919a' }; % R0425_20220728a debugging because the intan side was left on for 15 hours; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +target_Fs = 500; % in Hz, target LFP sampling rate after decimating the raw signal +convert_to_microvolts = false; + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; % changed probe_types.RatID to probe_types.ratID due to error + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); + session_dirs = dir(fullfile(rawdata_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + phys_folder = find_physiology_data(cur_dir); + + if isempty(phys_folder) + sprintf('no physiology data found for %s', session_name) + continue + end + + if any(strcmp(session_name, sessions_to_ignore)) % Jen added this in to ignore sessions as an attempt to debut "too many input arguments" + continue; + end + + lfp_fname = strcat(session_name, '_lfp.mat'); + processed_session_folder = fullfile(processed_folder, session_name); + full_lfp_name = fullfile(processed_session_folder, lfp_fname); + if ~isfolder(processed_session_folder) + mkdir(processed_session_folder) + end + if exist(lfp_fname) % the lfp_fname exists in the processed folder but this section is not recognizing it exists so still writes the file. + sprintf('%s already calculated, skipping', lfp_fname) + continue + end + + sprintf('working on %s', session_name) + [lfp, actual_Fs] = calculate_monopolar_LFPs_DL(phys_folder, target_Fs, convert_to_microvolts); + + save(full_lfp_name, 'lfp', 'actual_Fs', 'convert_to_microvolts'); + + end + +end \ No newline at end of file diff --git a/LFPs/script_extract_LFPs_DL.m b/LFPs/script_extract_LFPs_DL.m new file mode 100644 index 00000000..229cc8d6 --- /dev/null +++ b/LFPs/script_extract_LFPs_DL.m @@ -0,0 +1,72 @@ +% script to calculate LFPs for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a','R0427_20220919a' }; % R0425_20220728a debugging because the intan side was left on for 15 hours; + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +target_Fs = 500; % in Hz, target LFP sampling rate after decimating the raw signal +convert_to_microvolts = false; + +num_rats = length(ratIDs); + +for i_rat = 1 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; % changed probe_types.RatID to probe_types.ratID due to error + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); + session_dirs = dir(fullfile(rawdata_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + phys_folder = find_physiology_data(cur_dir); + + if isempty(phys_folder) + sprintf('no physiology data found for %s', session_name) + continue + end + + if any(strcmp(session_name, sessions_to_ignore)) % Jen added this in to ignore sessions as an attempt to debut "too many input arguments" + continue; + end + + lfp_fname = strcat(session_name, '_lfp.mat'); + processed_session_folder = fullfile(processed_folder, session_name); + full_lfp_name = fullfile(processed_session_folder, lfp_fname); + if ~isfolder(processed_session_folder) + mkdir(processed_session_folder) + end + if isfile(full_lfp_name) == 1 % the lfp_fname exists in the processed folder but this section is not recognizing it exists so still writes the file. + sprintf('%s already calculated, skipping', lfp_fname) + continue + end + + sprintf('working on %s', session_name) + [lfp, actual_Fs] = calculate_monopolar_LFPs_DL(phys_folder, target_Fs, convert_to_microvolts); + + save(full_lfp_name, 'lfp', 'actual_Fs', 'convert_to_microvolts'); + + end + +end \ No newline at end of file diff --git a/LFPs/script_extract_bipolar_LFPs_DL.m b/LFPs/script_extract_bipolar_LFPs_DL.m new file mode 100644 index 00000000..34490ea6 --- /dev/null +++ b/LFPs/script_extract_bipolar_LFPs_DL.m @@ -0,0 +1,63 @@ +% script to calculate LFPs for all of Jen's rats; store in files in +% the processed data folders + +parent_directory = 'X:\Neuro-Leventhal\data\ChoiceTask'; +summary_xls = 'ProbeSite_Mapping_MATLAB.xlsx'; +summary_xls_dir = 'X:\Neuro-Leventhal\data\ChoiceTask\Probe Histology Summary'; +summary_xls = fullfile(summary_xls_dir, summary_xls); + +probe_type_sheet = 'probe_type'; +probe_types = read_Jen_xls_summary(summary_xls, probe_type_sheet); +% NOTE - UPDATE FUNCTION read_Jen_xls_summary WHEN WE NEED OTHER +% INFORMATION OUT OF THAT SPREADSHEET + +[rat_nums, ratIDs, ratIDs_goodhisto] = get_rat_list(); + +target_Fs = 500; % in Hz, target LFP sampling rate after decimating the raw signal + +num_rats = length(ratIDs); + +for i_rat = 16 : num_rats + ratID = ratIDs{i_rat}; + rat_folder = fullfile(parent_directory, ratID); + + if ~isfolder(rat_folder) + continue; + end + + probe_type = probe_types{probe_types.ratID == ratID, 2}; + processed_folder = find_data_folder(ratID, 'processed', parent_directory); + rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); + session_dirs = dir(fullfile(processed_folder, strcat(ratID, '*'))); + num_sessions = length(session_dirs); + + for i_session = 1 : num_sessions + + session_name = session_dirs(i_session).name; + if strcmp(session_name, 'R0378_20210507a') || strcmp(session_name, 'R0425_20220728a') + % skip this session because for some reason only 63 channels + % were recorded + continue + end + cur_dir = fullfile(session_dirs(i_session).folder, session_name); + cd(cur_dir) + + lfp_fname = strcat(session_name, '_lfp.mat'); + full_lfp_name = fullfile(processed_folder, session_name, lfp_fname); + if isfile(full_lfp_name) + lfp_data = load(full_lfp_name); + else + continue + end + actual_Fs = lfp_data.actual_Fs; + + sprintf('working on %s', session_name) + lfp_bipolar_name = strcat(session_name, '_bipolar_lfp.mat'); + lfp_bipolar_name = fullfile(processed_folder, session_name, lfp_bipolar_name); + [bipolar_lfp, intan2probe_mapping] = calculate_bipolar_LFPs(lfp_data.lfp, probe_type); + + save(lfp_bipolar_name, 'bipolar_lfp', 'actual_Fs', 'probe_type', 'intan2probe_mapping', 'full_lfp_name'); + + end + +end \ No newline at end of file diff --git a/intan_choice_task_navigation_utils/create_scalo_folder.m b/intan_choice_task_navigation_utils/create_scalo_folder.m new file mode 100644 index 00000000..b7b9946e --- /dev/null +++ b/intan_choice_task_navigation_utils/create_scalo_folder.m @@ -0,0 +1,18 @@ +function scalo_folder = create_scalo_folder(session_name,event_name,parent_directory) +%UNTITLED2 Summary of this function goes here +% Detailed explanation goes here + +session_name_parts = split(session_name,'_'); +ratID = session_name_parts{1}; +processed_folder = find_data_folder(ratID, 'processed', parent_directory); + +session_folder = fullfile(processed_folder, session_name); + +if exist(session_folder, 'dir') + scalo_folder = sprintf('%s_scalos_%s', session_name, event_name); + scalo_folder = fullfile(session_folder, scalo_folder); +else + scalo_folder = ''; +end + +end \ No newline at end of file diff --git a/intan_choice_task_navigation_utils/find_data_folder.m b/intan_choice_task_navigation_utils/find_data_folder.m new file mode 100644 index 00000000..ca65c9e4 --- /dev/null +++ b/intan_choice_task_navigation_utils/find_data_folder.m @@ -0,0 +1,28 @@ +function session_path = find_data_folder(session_name, data_type, parent_directory) +%UNTITLED Summary of this function goes here +% Detailed explanation goes here + + +str_parts = split(session_name, '_'); +ratID = str_parts{1}; + +rd_folder = find_rat_data_folder(parent_directory, ratID, data_type); + +if length(ratID) == 5 + session_path = rd_folder; +elseif isfolder(fullfile(rd_folder, session_name)) + session_path = fullfile(rd_folder, session_name); +else + session_path = ''; +end + +end + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +function rd_folder = find_rat_data_folder(parent_directory, ratID, data_type) + +rat_folder = fullfile(parent_directory, ratID); +rd_folder = fullfile(rat_folder, strcat(ratID, '-', data_type)); + +end \ No newline at end of file diff --git a/intan_choice_task_navigation_utils/find_log_file.m b/intan_choice_task_navigation_utils/find_log_file.m new file mode 100644 index 00000000..6a55b203 --- /dev/null +++ b/intan_choice_task_navigation_utils/find_log_file.m @@ -0,0 +1,30 @@ +function log_fname = find_log_file(session_name, parent_directory) +%UNTITLED15 Summary of this function goes here +% Detailed explanation goes here + +session_name_parts = split(session_name, '_'); +ratID = session_name_parts{1}; +session_date = session_name_parts{2}(1:8); + +rawdata_folder = find_data_folder(ratID, 'rawdata', parent_directory); +rawdata_session_folder = fullfile(rawdata_folder, session_name); + +test_name = strcat(ratID, '_', session_date, '_*.log'); +test_name = fullfile(rawdata_session_folder, test_name); + +poss_logs = dir(test_name); + +log_fname = ''; +if isempty(poss_logs) + sprintf('no log file found for %s', session_name); +else + for ii = 1 : length(poss_logs) + if contains(poss_logs(ii).name, 'old') + continue + else + log_fname = fullfile(poss_logs(ii).folder, poss_logs(ii).name); + end + end +end + +end \ No newline at end of file diff --git a/intan_choice_task_navigation_utils/find_rawdata_folder.m b/intan_choice_task_navigation_utils/find_rawdata_folder.m new file mode 100644 index 00000000..dab30548 --- /dev/null +++ b/intan_choice_task_navigation_utils/find_rawdata_folder.m @@ -0,0 +1,67 @@ +function session_path = find_rawdata_folder(input1, varargin) +%UNTITLED Summary of this function goes here +% Detailed explanation goes here + +parent_directory = 'Z:\data\ChoiceTask'; + +if nargin > 1 + if strcmpi(varargin{1}, 'parent_directory') + parent_directory = varargin{2}; + end +end +if nargin > 2 + if strcmpi(varargin{2}, 'parent_directory') + parent_directory = varargin{3}; + end +end + +if nargin == 1 || (nargin > 1 && strcmpi(varargin{1}, 'parent_directory')) + % assume the full session name was given - i.e., 'RXXXX_YYYYMMDDz' + + str_parts = split(input1, '_'); + ratID = str_parts{1}; + + rat_rd_folder = find_rat_rawdata_folder(parent_directory, ratID); + + if isfolder(fullfile(rat_rd_folder, input1)) + session_path = fullfile(rat_rd_folder, input1); + else + session_path = ''; + end + +elseif nargin == 2 + % assume the first argument is the rat number or ratID, the second + % argument is the date string or date (if date, assume it's session a) + + if isinteger(input1) + ratID = sprintf('R%04d', input1); + elseif ischar(input1) + ratID = input1; + end + rat_rd_folder = find_rat_rawdata_folder(parent_directory, ratID); + + if ischar(varargin{1}) + session_name = strcat(ratID, '_', varargin{1}); + elseif isdatetime(varargin{1}) + date_string = datestr(varargin{1}, 'yyyymmdd'); + session_name = strcat(ratID, '_', date_string, 'a'); + end + + if isfolder(fullfile(rat_rd_folder, session_name)) + session_path = fullfile(rat_rd_folder, session_name); + else + session_path = ''; + end + +end + +end + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +function rat_rd_folder = find_rat_rawdata_folder(parent_directory, ratID) + +rat_folder = fullfile(parent_directory, ratID); +rat_rd_folder = fullfile(rat_folder, strcat(ratID, '-rawdata')); + +end \ No newline at end of file diff --git a/intan_choice_task_navigation_utils/get_rawdata_ephys_folder.m b/intan_choice_task_navigation_utils/get_rawdata_ephys_folder.m new file mode 100644 index 00000000..6f23d2b4 --- /dev/null +++ b/intan_choice_task_navigation_utils/get_rawdata_ephys_folder.m @@ -0,0 +1,26 @@ +function rawdata_ephys_folder = get_rawdata_ephys_folder(rawdata_folder,session_name) +%UNTITLED13 Summary of this function goes here +% Detailed explanation goes here + +rawdata_session_folder = fullfile(rawdata_folder, session_name); + +session_name_parts = split(session_name, '_'); +session_date_str = session_name_parts{2}(1:8); +intan_date_str = session_date_str(3:8); + +test_name = strcat(session_name_parts{1}, '_', session_date_str, '*', intan_date_str, '*'); +test_name = fullfile(rawdata_session_folder, test_name); + +rawdata_ephys_dirs = dir(test_name); + +if length(rawdata_ephys_dirs) == 1 + rawdata_ephys_folder = fullfile(rawdata_ephys_dirs.folder, rawdata_ephys_dirs.name); +elseif length(rawdata_ephys_dirs) > 1 + sprintf('multiple ephys folders found for %s', session_name) + rawdata_ephys_folder = ''; +else + sprintf('no ephys folders found for %s', session_name) + rawdata_ephys_folder = ''; +end + +end \ No newline at end of file diff --git a/silicon_probe_LFP_analysis/calculate_monopolar_LFPs_DL.m b/silicon_probe_LFP_analysis/calculate_monopolar_LFPs_DL.m new file mode 100644 index 00000000..be3054b4 --- /dev/null +++ b/silicon_probe_LFP_analysis/calculate_monopolar_LFPs_DL.m @@ -0,0 +1,113 @@ +function [lfp_data, actual_lfpFs] = calculate_monopolar_LFPs_DL(intan_folder, target_Fs, convert_to_microvolts) +% +% INPUTS +% intan_folder - folder containing intan data (amplifier.dat, info.rhd) +% target_Fs - target sampling rate + +raw_block_size = 100000; % number of samples to handle at a time (titrate to memory), may want to make this a varargin +bytes_per_sample = 2; +filtOrder = 1000; + +cd(intan_folder); + +rhd_file = dir('*.rhd'); +if length(rhd_file) > 1 + error('more than one rhd file in ' + intan_folder); + return +elseif isempty(rhd_file) + error('no rhd files found in ' + intan_folder); + return +end + +amp_file = dir('amplifier.dat'); %rename for Watson Data; change back to amplifier.dat +if isempty(amp_file) + error('no amplifier files found in ' + intan_folder); + return +end + +rhd_info = read_Intan_RHD2000_file_DL(rhd_file.name); +amplifier_channels = rhd_info.amplifier_channels; +num_channels = length(amplifier_channels); +Fs = rhd_info.frequency_parameters.amplifier_sample_rate; + +samples_per_channel = amp_file.bytes / (num_channels * bytes_per_sample); + +r = round(Fs / target_Fs); +actual_lfpFs = Fs/r; +raw_overlap_size = ceil(filtOrder * 2 / r) * r; % can use this line to check the overlap by making it equal to zero +lfp_overlap_size = raw_overlap_size / r; + +lfp_block_size = raw_block_size / r; +num_lfp_samples = ceil(samples_per_channel / r); +lfp_data = zeros(num_channels, num_lfp_samples); + +% calculate the number of blocks that will be needed +net_lfp_samples_per_block = lfp_block_size - lfp_overlap_size; +num_blocks = ceil(num_lfp_samples / net_lfp_samples_per_block); + +currentLFP = zeros(num_channels, lfp_block_size); + +% do the first block separate from the rest, since there will be some +% overlap for the rest of the blocks +LFPstart = 1; +LFPend = (LFPstart + lfp_block_size - 1) - lfp_overlap_size; + +disp(['Block 1 of ' num2str(num_blocks)]); +amplifier_data = readIntanAmplifierData_by_sample_number(amp_file.name,1,raw_block_size,amplifier_channels,convert_to_microvolts); +for i_ch = 1 : num_channels + currentLFP(i_ch, :) = ... + decimate(amplifier_data(i_ch, :), r, filtOrder, 'fir'); +end +lfp_data(:, LFPstart:LFPend) = currentLFP(:, 1:LFPend); +clear currentLFP; + +LFPstart = LFPend + 1; +raw_block_plus_overlap_size = raw_block_size + raw_overlap_size; +LFPblock_start = lfp_overlap_size + 1; +LFPblock_end = (LFPblock_start + lfp_block_size - 1) - lfp_overlap_size; + +read_final_samples = false; +for i_block = 2 : num_blocks + + disp(['Block ' num2str(i_block) ' of ' num2str(num_blocks)]); + +% read_start_sample = (i_block-1) * raw_block_size - raw_overlap_size + 1; + read_start_sample = (LFPstart-1) * r - raw_overlap_size; + read_end_sample = read_start_sample + raw_block_plus_overlap_size - 1; + + new_amplifier_data = readIntanAmplifierData_by_sample_number(amp_file.name,read_start_sample,read_end_sample,amplifier_channels,convert_to_microvolts); + + if i_block < num_blocks + if read_end_sample > samples_per_channel + % rarely, the end of the padded block goes past the end of the + % file, and messes up the indices + clear currentLFP + LFPend = size(lfp_data, 2); + read_final_samples = true; + else + LFPend = (LFPstart + lfp_block_size - 1) - lfp_overlap_size; + end + else + clear currentLFP; + read_final_samples = true; + LFPend = size(lfp_data, 2); + end + + for i_ch = 1 : num_channels + currentLFP(i_ch,:) = ... + decimate(new_amplifier_data(i_ch, :), r, filtOrder, 'fir'); + end + + if read_final_samples + LFPblock_end = size(currentLFP, 2); + end + try + lfp_data(:, LFPstart:LFPend) = currentLFP(:, LFPblock_start : LFPblock_end); + catch + keyboard + end + if read_final_samples + break + end + LFPstart = LFPend + 1; +end \ No newline at end of file diff --git a/silicon_probe_LFP_analysis/script_extract_power_spectra.m b/silicon_probe_LFP_analysis/script_extract_power_spectra.m index 81af1faa..8202348b 100644 --- a/silicon_probe_LFP_analysis/script_extract_power_spectra.m +++ b/silicon_probe_LFP_analysis/script_extract_power_spectra.m @@ -14,9 +14,11 @@ NN8x8 = ["R0326", "R0327", "R0372", "R0379", "R0374", "R0376", "R0378", "R0394", "R0395", "R0396", "R0412", "R0413"]; % Specify list of ratID associated with each probe_type ASSY156 = ["R0411", "R0419"]; -ASSY236 = ["R0420", "R0425", "R0427", "R0457"]; +ASSY236 = ["R0420", "R0425", "R0427", "R0457", "R0456", "R0459", "R0460", "R0463", "R0465", "R0466"]; -sessions_to_ignore = {'R0378_20210507a', 'R0425_20220728a', 'R0427_20220920a'}; % R0425_20220728a debugging because the intan side was left on for 15 hours; +sessions_to_ignore = {'R0378_20210507a', 'R0326_20191107a', 'R0425_20220728a', 'R0425_20220816b', 'R0427_20220920a', 'R0427_20220920a'}; +sessions_to_ignore1 = {'R0425_20220728_ChVE_220728_112601', 'R0427_20220920_Testing_220920_150255'}; +sessions_to_ignore2 = {'R0427_20220908a', 'R0427_20220909a', 'R0427_20220912a','R0427_20220913a', 'R0427_20220914a', 'R0427_20220915a', 'R0427_20220916a'}; % R0427_20220920a does not have an 'info.rhd' file %% @@ -56,19 +58,46 @@ lfp_file = dir(fullfile(session_path_processed, '**', '*_lfp.mat')); cd(session_path_processed); + + for i_file = 1 : length(lfp_file) + lfp_fname = lfp_file(i_file).name; + name_parts = split(lfp_fname, '_'); + + if strcmp(name_parts(3), 'bipolar') + bipolar_lfpname = fullfile(lfp_file(i_file).folder, lfp_fname); + else + monopolar_lfpname = fullfile(lfp_file(i_file).folder, lfp_fname); + end + end + % lfp = load(lfp_fname.name); % I think the lfp needs to be loaded % in the lfp_NNsite_order script in the next line of this fxn. % LFP file needs to be loaded before the [power_lfps,f] function - lfp_fname = fullfile(lfp_file.folder, lfp_file.name); +% lfp_fname = fullfile(lfp_file.folder, lfp_file.name); -% if exist(power_fn, 'file') && exist(diff_power_fn, 'file') -% continue -% end + if exist(power_fn, 'file') && exist(diff_power_fn, 'file') + continue + end - lfp_data = load(lfp_fname); - Fs = lfp_data.actual_Fs; + for i_lfp = 1 : 2 + if i_lfp == 1 + isbipolar = true; + else + isbipolar = false; + end + + if isbipolar + lfp_fname = bipolar_lfpname; + else + lfp_fname = monopolar_lfpname; + end + + lfp_data = load(lfp_fname); + Fs = lfp_data.actual_Fs; + end + num_rows = size(lfp_data.lfp,1); % for now, skipping R0378_20210507a because the session only recorded 63 channels instead of 64. Need to rewrite lfp_NNsite_order and diff functions to fix this issue by determining which channel was not recorded. if num_rows < 64 continue