123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366 |
- %%
- clear
- general_directory = '/Users/b1037210/Desktop/Progetti/JA_Stage1/'
- cd(general_directory)
- addpath(general_directory)
- %%
- addpath(genpath('/Users/b1037210/Documents/MATLAB/fieldtrip-20161117'))
- ft_defaults
- %%
- subjects = ...
- {'S01_pilot',...
- 'S02_pilot',...
- 'S03_pilot',...
- 'S04_pilot',...
- 'S05_pilot',...
- 'S06_pilot',...
- 'S07_pilot',...
- 'S08_pilot',...
- 'S09_pilot',...
- 'S10_pilot',...
- 'S11_pilot',...
- 'S12_pilot'}
- threshold = {[0.018 0.032],...
- [0.016 0.035],...
- [0.018 0.040],...
- [0.019 0.039],...
- [0.017 0.037],...
- [0.016 0.036],...
- [0.019 0.042],...
- [0.018 0.039],...
- [0.016 0.035],...
- [0.015 0.035],...
- [0.02 0.03],...
- [0.02 0.04]}
- TMS_timings = {'-0.05', '0.0'};
- final_table = table;
- %%
- for s = 1:length(subjects)
-
- %% cd directory
- cd([general_directory,subjects{s}])
-
- %% ...for each tms timing
- for t = 1:length(TMS_timings)
-
- %% file loading
- file_tmp = (dir(strcat(subjects{s},'_', TMS_timings{t},'_TMS*')))
-
-
- if size(file_tmp,1)>1
-
- disp(strcat('concatenate ',file_tmp(1).name))
- [y SR MetaData] = GB_gtec_files_concatenation(file_tmp);
-
- else
-
- disp(strcat('loading ', file_tmp.name))
- load(file_tmp.name)
- end
-
- clear file_tmp
- %% name the trigger channel and time channel, then erase time from y
-
- % select trigger row and time row
- trigger = y(end,:); % trigger channel for glove touch detection
- time = y(1,:); % time-series in seconds
-
- % new data without time, so channels in the gtect settings and data
- % analysis have the same nnumber
- y = y(2:end,:);
-
- %% get triggers and sample points at which each event of interest happens
-
- % get were variations happen within trigger channel
- diff_trigger = diff(trigger);
-
- % get trigger codes wihtin trgger line
- trigger_numbers = unique(trigger);
-
- % find were the wash_out ends
- end_wash_out = find(diff_trigger == 1)+1;
- end_wash_out = end_wash_out(1);
-
- % for each code present, get the sample were the rising triggers are present
- for x = 1:length(trigger_numbers)
-
- % find rising triggers and add 1 sampling point (due to diff which removes one)
- tmp = find(diff_trigger == trigger_numbers(x))+1;
-
- % based on the code --> assign the tmp numbers to the correct variable
- switch trigger_numbers(x)
-
- % case 14 % color press CF
- % sp_slow_id = tmp(tmp>end_wash_out);
- %
- % case 15 % color press CF
- % sp_fast_id = tmp(tmp>end_wash_out);
-
- case 16 % color press CF
- sp_out_slow_id = tmp(tmp>end_wash_out);
-
- case 17 % color press CF
- sp_out_fast_id = tmp(tmp>end_wash_out);
-
- % case 18 % color press CF
- % sp_out_ok_catch_id = tmp(tmp>end_wash_out);
- %
- % case 19 % color press CF
- % sp_out_failed_catch_id = tmp(tmp>end_wash_out);
- %
-
-
-
-
- case 20 % color press CF
- sp_color_CF_press = tmp(tmp>end_wash_out);
-
- case 21 % color lift CF
- sp_color_CF_lift = tmp(tmp>end_wash_out);
-
- case 22 % color catch CF
- sp_color_CF_catch = tmp(tmp>end_wash_out);
-
- case 30 % color press SBJ
- sp_color_SBJ_press = tmp(tmp>end_wash_out);
-
- case 31 % color_lift SBJ
- sp_color_SBJ_lift = tmp(tmp>end_wash_out);
-
- case 40 % CF hitting
- sp_CF_hitting = tmp(tmp>end_wash_out);
-
- case 41 % SBJ hitting
- sp_SBJ_hitting = tmp(tmp>end_wash_out);
-
- case 128 % TMS pulse
- sp_TMS = tmp(tmp>end_wash_out);
- end
-
- end
-
- %%
-
- % check dimensions
- size(sp_color_CF_press)
- size(sp_color_CF_lift)
- size(sp_color_CF_catch)
- size(sp_TMS)
-
- %
- [CF_colors_idx] = sortrows(vertcat(sp_color_CF_press',sp_color_CF_lift',sp_color_CF_catch'));
- CF_colors_order = CF_colors_idx;
-
- % assign codes to sample points corresponding to press lift and catch
- CF_colors_order(ismember(CF_colors_idx, sp_color_CF_press)) = 1;
- CF_colors_order(ismember(CF_colors_idx, sp_color_CF_lift)) = -1;
- CF_colors_order(ismember(CF_colors_idx, sp_color_CF_catch)) = 0;
-
- %%
-
- size(sp_out_slow_id)
- size(sp_out_fast_id)
-
- [trials_idx] = sortrows(vertcat(sp_out_slow_id',sp_out_fast_id'))
-
-
- %%
- for ex = 1:length(trials_idx)
-
- [v_abs(ex) p_abs(ex)] = min(abs(CF_colors_idx - trials_idx(ex)))
- ddd = CF_colors_idx - trials_idx(ex);
- exclude(ex) = p_abs(ex);
- is_exclude_positive(ex) = ddd(exclude(ex))
- end
-
- clear v_abs
- clear ddd
- clear is_exclude_positive
-
-
-
- %% EMG preprocessing
-
- ECD_right = y(2,:);
- FDS_right = y(4,:);
- ECD_response = y(6,:);
-
-
-
- Fsample = SR;
- hp_filter = 10
- lp_filter = 2500
-
-
- ECD_right_filt = ft_preproc_highpassfilter(ECD_right, Fsample, hp_filter,1);
- FDS_right_filt = ft_preproc_highpassfilter(FDS_right, Fsample, hp_filter,1);
-
- ECD_right_filt = ft_preproc_lowpassfilter(ECD_right_filt, Fsample, lp_filter,1);
- FDS_right_filt = ft_preproc_lowpassfilter(FDS_right_filt, Fsample, lp_filter,1);
-
-
- ECD_right_filt = ft_preproc_dftfilter(ECD_right_filt, Fsample);
- FDS_right_filt = ft_preproc_dftfilter(FDS_right_filt, Fsample);
-
-
- %% take away the decay from MEPs
-
- cfg = [];
- cfg.start_time = 0.005;
- cfg.end_time = 0.15;
- cfg.SR = SR;
- cfg.sp_trigger = sp_TMS;
- [ECD_sweep_cuts_tmp, ECD_sweep_time, ECD_sample_points] = gb_cutting_sweeps(cfg, ECD_right_filt);
- [FDS_sweep_cuts_tmp, FDS_sweep_time, FDS_sample_points] = gb_cutting_sweeps(cfg, FDS_right_filt);
- %% get MEPs p2p amplitudes
- cfg = [];
- cfg.start_time = threshold{s}(1);
- cfg.end_time = threshold{s}(2);
- cfg.SR = SR;
- cfg.sp_trigger = sp_TMS;
-
- [ECD_MEP_amplitudes] = gb_MEP_p2p(cfg, ECD_right_filt);
- [FDS_MEP_amplitudes] = gb_MEP_p2p(cfg, FDS_right_filt);
-
-
- %% fit the decay and remove it
- cfg = [];
- cfg.start_time = -0.15;
- cfg.end_time = -0.002;
- cfg.SR = SR;
- cfg.sp_trigger = sp_TMS;
- [ECD_sweep_cuts_bkgrnd_tmp, ECD_sweep_time_bkgrnd, ECD_sample_points_bkgrnd] = gb_cutting_sweeps(cfg, ECD_right_filt)
- [FDS_sweep_cuts_bkgrnd_tmp, FDS_sweep_time_bkgrnd, FDS_sample_points_bkgrnd] = gb_cutting_sweeps(cfg, FDS_right_filt)
- %%
- for x = 1:length(ECD_sweep_cuts_bkgrnd_tmp)
- [~, ~, resECD{x}] = fit(ECD_sweep_time_bkgrnd',ECD_sweep_cuts_bkgrnd_tmp{x}','poly3'); % time in ms
- ECD_sweep_cuts_bkgrnd{x} =resECD{x}.residuals';
- [~, ~, resFDS{x}] = fit(FDS_sweep_time_bkgrnd',FDS_sweep_cuts_bkgrnd_tmp{x}','poly3'); % time in ms
- FDS_sweep_cuts_bkgrnd{x} =resFDS{x}.residuals';
- end
- clear resECD resFDS
- %%
- for x = 1:length(ECD_sample_points_bkgrnd)
-
- ECD_right_filt(ECD_sample_points_bkgrnd{x}) = ECD_sweep_cuts_bkgrnd{x};
- FDS_right_filt(FDS_sample_points_bkgrnd{x}) = FDS_sweep_cuts_bkgrnd{x};
- end
- %% any contraction befor TMS?
-
- cfg = [];
- cfg.start_time = -0.05;
- cfg.end_time = -0.002;
- cfg.SR = SR;
- cfg.sp_trigger = sp_TMS;
- cfg.analysis = 'p2p';
-
- [ECD_MEP_bkgrnd] = gb_background_activity(cfg, ECD_right_filt)
- [FDS_MEP_bkgrnd] = gb_background_activity(cfg, FDS_right_filt)
-
-
- %% muscle contraction before TMS?
-
- MEP_excl_ECD_idx = find(ECD_MEP_bkgrnd>200) % microvolt
- MEP_excl_FDS_idx = find(FDS_MEP_bkgrnd>200) % microvolt
-
-
- %% EXCLUDE MEP<50 p2p microVolts
-
- MEP_excl2_ECD_idx = find(ECD_MEP_amplitudes<50) % microvolt
-
-
- %% EXCLUDE out slow and out fast trials
-
- ECD_MEP_amplitudes(exclude) = nan;
-
-
- %% exclude trials with background contraction
-
- ECD_MEP_amplitudes(union(MEP_excl_ECD_idx, MEP_excl_FDS_idx)) = nan;
-
-
- ECD_MEP_amplitudes(MEP_excl2_ECD_idx) = nan;
-
-
- %% median amplitude based on color trigger
- summary.ECD_press{s,t} = nanmedian((ECD_MEP_amplitudes(CF_colors_order == 1)));
- summary.ECD_lift{s,t} = nanmedian((ECD_MEP_amplitudes(CF_colors_order == -1)));
- %%
-
- cfg = [];
- cfg.start_time = -0.015;
- cfg.end_time = 0.15;
- cfg.SR = SR;
- cfg.sp_trigger = sp_TMS;
- [bbb.ECD_sweep_PLOT_CHECKS, bbb.ECD_sweep_time_PLOT_CHECKS, bbb.ECD_sample_points_PLOT_CHECKS] = gb_cutting_sweeps(cfg, ECD_right_filt);
-
- %%
-
- trial = (1:length(CF_colors_order))';
- subject_rep = repmat(subjects(s), [length(CF_colors_order),1]);
- TMS_timings_rep = repmat(TMS_timings(t), [length(CF_colors_order),1]);
-
- subject_table{t} = table(subject_rep,trial, ECD_MEP_amplitudes,...
- ECD_MEP_bkgrnd,...
- FDS_MEP_bkgrnd,...
- TMS_timings_rep,...
- CF_colors_order)
-
-
-
-
-
-
- clearvars -except s t subjects TMS_timings threshold subject_table general_directory final_table summary
-
- %%
-
- end
-
- final_subject_table = vertcat(subject_table{:})
-
-
- output_filename = strcat('final_subject_table_',subjects{s},'_' , '.mat');
- save(output_filename,'final_subject_table')
-
-
- final_table = vertcat(final_table, final_subject_table)
-
- clearvars -except s t subjects TMS_timings threshold general_directory final_table summary
-
-
- end
- %%
- cd(general_directory)
- filename = 'final_table_for_R_10-2500_def_o1_rep.xlsx';
- writetable(final_table,filename,'Sheet',1)
|