Explorar el Código

Dateien hochladen nach ''

Johannes Wetekam hace 10 meses
padre
commit
41efdbe5e0
Se han modificado 1 ficheros con 729 adiciones y 0 borrados
  1. 729 0
      DD_Plotting.m

+ 729 - 0
DD_Plotting.m

@@ -0,0 +1,729 @@
+%% define some variables
+
+filt = 2; % 1: non, 2: 0.1 Hz, 3: 10 Hz, 4: 150 Hz, 5: 300 Hz
+cutTime = 1; % plot only certain time window? (0: no, 1: short, 2: long)
+avType = 1; % plot grand average (1) or individual averages (2)
+fixedY = 1; % use dynamic (0) or fixed (1) y-axis?
+inclCtrl = 2; % plot control data? (0: no, 1: 50 %, 2: MR, 3: both controls)
+statsType = 2; % decide which statistics to use if control is included; 1: t-test, 2: ANOVA
+plotWinI = 1; % plot RMS-windows?
+plotRespI = 1; % plot responses?
+plotDifI = 0; % plot difference curves?
+plotMtrxI = 0; % plot all responses in a matrix?
+plotHmI = 0; % plot heatmap?
+hmDataI = 1; % choose data that should be plotted as heatmap (1: Cohens D; 2: RMS of control, only useful for 50Per dataset)
+plotBpI = 0; % plot boxplots?
+stimMark = 0; % decide how to plot the stimulus (no: 0, bar: 1)
+cStart = 1; % stimulus combination to start with
+oDisEloc = 0; % plot only DisNoAM and Eloc?
+saveI = 0; % save figure?
+
+% decide which set of rms-windows to plot
+switch filt
+        case {1,2,3}
+            switch cutTime
+                case 0
+                    uWin = (1:5);
+                    boxpFigWi = 800; % width of boxplot figure
+                case 1
+                    uWin = (4:4);
+                    boxpFigWi = 300; % width of boxplot figure
+                case 2
+                    uWin = (4:5);
+                    boxpFigWi = 500; % width of boxplot figure
+            end
+    case {4,5}
+        uWin = (1:3);
+        boxpFigWi = 500; % width of boxplot figure
+end
+% uWin = uWin+5;
+
+switch filt
+    case 1
+        filtName = 'noF';
+        y_lim = [-3.5,13];
+    case 2
+        filtName = '0.1F';
+        y_lim = [-3.5,13];
+    case 3
+        filtName = '10F';
+        y_lim = [-3.5,13];
+    case 4
+        filtName = '150F';
+        y_lim = [-4.5,4.5];
+    case 5
+        filtName = '300F';
+        y_lim = [-4,4];
+end
+
+% define colors and names
+colDev = [1,0,0]; % dev: red
+colStd = [0,0,1]; % std: blue
+col50Per = [1,0,1]; % 50 %: magenta
+colMR = [0,0,0]; % MR: black
+nameDev = 'Deviant';
+nameStd = 'Standard';
+name50Per = '50 % Control';
+nameMR = 'MS Control';
+
+%% load processed data
+
+switch inclCtrl
+    case 0
+        load('DD_statsData_noCtrl.mat') % load stats data
+    case 1 % load 50Per data and rename important variables
+        load('DD_avData_50Per.mat')
+        dataAv_cell_Ctrl = dataAv_cell; % rename
+        dataGrAv_cell_Ctrl = dataGrAv_cell;
+        dataDifAv_cell_Ctrl = dataDifAv_cell; % rename
+        dataDifGrAv_cell_Ctrl = dataDifGrAv_cell;
+        dataSe_cell_Ctrl = dataSe_cell;
+        filenames_cell_Ctrl = filenames_cell;
+        load('DD_statsData_50Per.mat') % load stats data
+        nameCtrl = name50Per;
+        colCtrl = col50Per; % set control color to 50 % color
+    case 2 % load MR data and rename important variables
+        load('DD_avData_MR.mat')
+        dataAv_cell_Ctrl = dataAv_cell; % rename
+        dataGrAv_cell_Ctrl = dataGrAv_cell;
+        dataDifAv_cell_Ctrl = dataDifAv_cell; % rename
+        dataDifGrAv_cell_Ctrl = dataDifGrAv_cell;
+        dataSe_cell_Ctrl = dataSe_cell;
+        filenames_cell_Ctrl = filenames_cell;
+        load('DD_statsData_MR.mat') % load stats data
+        nameCtrl = nameMR;
+        colCtrl = colMR; % set control color to MR color
+    case 3
+        load('DD_avData_50Per.mat')
+        dataAv_cell_50Per = dataAv_cell; % rename
+        dataGrAv_cell_50Per = dataGrAv_cell;
+        dataDifAv_cell_50Per = dataDifAv_cell; % rename
+        dataDifGrAv_cell_50Per = dataDifGrAv_cell;
+        dataSe_cell_50Per = dataSe_cell;
+        filenames_cell_50Per = filenames_cell;
+        load('DD_avData_MR.mat')
+        dataAv_cell_MR = dataAv_cell; % rename
+        dataGrAv_cell_MR = dataGrAv_cell;
+        dataDifAv_cell_MR = dataDifAv_cell; % rename
+        dataDifGrAv_cell_MR = dataDifGrAv_cell;
+        dataSe_cell_MR = dataSe_cell;
+        filenames_cell_MR = filenames_cell;
+        load('DD_statsData_bothCtrl.mat') % load stats data
+end
+
+% load Oddball data (will overwrite Ctrl variables that were not renamed)
+load('DD_avData_Oddball.mat')
+%% do some additional calculations
+
+switch zsNormI
+    case 0
+        ylableName = 'Voltage [µV]';
+    case 1
+        ylableName = 'z-norm. voltage [a.u.]';
+end
+
+switch cutTime
+    case 1
+        switch filt
+            case {1,2,3}
+                timeStart = -5; % starting time of the plot (relative to stim onset)
+                timeEnd = 15; % ending time of the plot (relative to stim onset)
+            case {4,5}
+                timeStart = -5; % starting time of the plot (relative to stim onset)
+                timeEnd = 15; % ending time of the plot (relative to stim onset)
+        end
+    case 2
+        timeStart = -5; % starting time of the plot (relative to stim onset)
+        timeEnd = 40; % ending time of the plot (relative to stim onset)        
+end
+
+% some calculations for matrix- and heatmap-plotting
+combNameS = split(combName(cStart:end),","); % split comb names at comma
+if plotMtrxI==1||plotHmI==1
+    oDisEloc = 0;
+    cStart = 2; % skip first stimulus combination (pure tones)
+    stimOrd = ["DisNoAM","DisAM","DisMimic","Eloc","ElocMimic"]; % define order of stimuli to plot in the matrix
+    nStim = size(stimOrd,2);
+    combCoord = zeros(nComb-1,2); % preallocate
+    for c = 1:nComb-1
+        for s = 1:2
+            combCoord(c,s) = find(ismember(stimOrd,combNameS(c,s)));
+        end
+    end
+    nTiles = 1:nStim*nStim; % number of tiles in matrix
+    tileMtrx = reshape(nTiles,[nStim,nStim])'; % create matrix to use for indexing later
+end
+
+if oDisEloc==1 % plot only stim combination 7 (DisNoAM/Eloc) if oDisEloc is 1
+    cStart = 7;
+    nComb = 7;
+end
+%% plotting: responses
+
+if plotRespI==1
+    switch plotMtrxI
+        case 1
+            figure('NumberTitle','off','Name','All Stims','Position',[0,0,2000,1300],'Renderer','painters')
+            t = tiledlayout(nStim,nStim);
+    end
+    for c = cStart:nComb % run once for each stimulus combination
+
+        % extract data corresponding to current stimulation-condition from
+        % cell-arrays
+        switch plotDifI
+            case 0
+                dataAv = dataAv_cell{c};
+                dataGrAv = dataGrAv_cell{c};
+                dataSe = dataSe_cell{c};
+                switch inclCtrl
+                    case {1,2}
+                        dataAvCtrl = dataAv_cell_Ctrl{c};
+                        dataGrAvCtrl = dataGrAv_cell_Ctrl{c};
+                        dataSeCtrl = dataSe_cell_Ctrl{c};
+                        filenames_Ctrl = filenames_cell_Ctrl{c};
+                        bsPos = cell2mat(strfind(filenames_Ctrl,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_Ctrl = extractAfter(filenames_Ctrl,startPos);
+                    case 3
+                        dataAv50Per = dataAv_cell_50Per{c};
+                        dataGrAv50Per = dataGrAv_cell_50Per{c};
+                        dataSe50Per = dataSe_cell_50Per{c};
+                        filenames_50Per = filenames_cell_50Per{c};
+                        bsPos = cell2mat(strfind(filenames_50Per,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_50Per = extractAfter(filenames_50Per,startPos);
+                        dataAvMR = dataAv_cell_MR{c};
+                        dataGrAvMR = dataGrAv_cell_MR{c};
+                        dataSeMR = dataSe_cell_MR{c};
+                        filenames_MR = filenames_cell_MR{c};
+                        bsPos = cell2mat(strfind(filenames_MR,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_MR = extractAfter(filenames_MR,startPos);
+                end
+            case 1
+                dataAv = dataDifAv_cell_Ctrl{c};
+                dataGrAv = dataDifGrAv_cell_Ctrl{c};
+                plotWinI = 0; % plot no significance windows when plotting diference curves
+                switch inclCtrl
+                    case 0
+                        error('To plot difference curves, activate inclCtrl')
+                    case {1,2}
+                        dataAvCtrl = dataAv_cell_Ctrl{c}; % won't be plotted, only processed to keep the code simpler
+                        dataGrAvCtrl = dataGrAv_cell_Ctrl{c};
+                        filenames_Ctrl = filenames_cell_Ctrl{c};
+                        bsPos = cell2mat(strfind(filenames_Ctrl,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_Ctrl = extractAfter(filenames_Ctrl,startPos);
+                    case 3
+                        dataAv50Per = dataAv_cell_50Per{c}; % won't be plotted, only processed to keep the code simpler
+                        dataGrAv50Per = dataGrAv_cell_50Per{c};
+                        filenames_50Per = filenames_cell_50Per{c};
+                        bsPos = cell2mat(strfind(filenames_50Per,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_50Per = extractAfter(filenames_50Per,startPos);
+                        dataAvMR = dataAv_cell_MR{c}; % won't be plotted, only processed to keep the code simpler
+                        dataGrAvMR = dataGrAv_cell_MR{c};
+                        filenames_MR = filenames_cell_MR{c};
+                        bsPos = cell2mat(strfind(filenames_MR,'\')); % find positions of back slashes in filename
+                        startPos = bsPos(:,end); % use last back slash position as starting point
+                        filenames_MR = extractAfter(filenames_MR,startPos);
+                end
+        end
+        filenames = filenames_cell{c};
+        recID = recID_cell{c};
+        bsPos = cell2mat(strfind(filenames,'\')); % find positions of back slashes in filename
+        startPos = bsPos(:,end); % use last back slash position as starting point
+        filenames = extractAfter(filenames,startPos);
+        timetr = round(timetr_cell{c}*1000,4);
+        stimDur = stimDur_cell{c}(1)*1000;
+        stimDelay = stimDelay_cell{c}(1)*1000;
+        stimWin = [stimDelay,stimDelay+stimDur];
+        nFiles = nFiles_cell{c};
+        fsDwn = fsDwn_cell{c};
+
+        % process timetrace
+        switch cutTime
+            case 0
+                timeWin = floor([timetr(1),timetr(end)]); % data in this window will be plotted (relative to recording onset)
+                timeStart = timeWin(1)-(stimDelay);
+                timeEnd = timeWin(2)-(stimDelay);
+                timeCut = (1:size(timetr,2));
+            case {1,2}
+                timeWin = [round((stimDelay)+timeStart,4),...
+                    round((stimDelay)+timeEnd,4)]; % data in this window will be plotted (relative to recording onset)
+                timeCut = round(timeWin(1)/1000*fsDwn:timeWin(2)/1000*fsDwn);
+        end
+
+        % subdivide responses for simplicity
+        ADevAv = dataAv{1};
+        AStdAv = dataAv{2};
+        BDevAv = dataAv{3};
+        BStdAv = dataAv{4};
+        ADevGrAv = dataGrAv{1};
+        AStdGrAv = dataGrAv{2};
+        BDevGrAv = dataGrAv{3};
+        BStdGrAv = dataGrAv{4};
+        ADevSe = dataSe{1};
+        AStdSe = dataSe{2};
+        BDevSe = dataSe{3};
+        BStdSe = dataSe{4};
+        switch inclCtrl
+            case {1,2}
+                AAvCtrl = dataAvCtrl{1};
+                BAvCtrl = dataAvCtrl{2};
+                AGrAvCtrl = dataGrAvCtrl{1};
+                BGrAvCtrl = dataGrAvCtrl{2};
+                ASeCtrl = dataSeCtrl{1};
+                BSeCtrl = dataSeCtrl{2};
+            case 3
+                AAv50Per = dataAv50Per{1};
+                BAv50Per = dataAv50Per{2};
+                AGrAv50Per = dataGrAv50Per{1};
+                BGrAv50Per = dataGrAv50Per{2};
+                ASe50Per = dataSe50Per{1};
+                BSe50Per = dataSe50Per{2};
+                AAvMR = dataAvMR{1};
+                BAvMR = dataAvMR{2};
+                AGrAvMR = dataGrAvMR{1};
+                BGrAvMR = dataGrAvMR{2};
+                ASeMR = dataSeMR{1};
+                BSeMR = dataSeMR{2};
+        end
+
+
+        % plot single channel
+        for s = 1:2 % once for A, once for B
+            switch s
+                case 1
+                    grAvStd = AStdGrAv;
+                    grAvDev = ADevGrAv;
+                    avStd = AStdAv;
+                    avDev = ADevAv;
+                    seStd = AStdSe;
+                    seDev = ADevSe;
+                    stimName = 'A';
+                    switch inclCtrl
+                        case {1,2}
+                            grAvCtrl = AGrAvCtrl;
+                            avCtrl = AAvCtrl;
+                            seCtrl = ASeCtrl;
+                        case 3
+                            grAv50Per = AGrAv50Per;
+                            av50Per = AAv50Per;
+                            se50Per = ASe50Per;
+                            grAvMR = AGrAvMR;
+                            avMR = AAvMR;
+                            seMR = ASeMR;
+                    end
+                case 2
+                    grAvStd = BStdGrAv;
+                    grAvDev = BDevGrAv;
+                    avStd = BStdAv;
+                    avDev = BDevAv;
+                    seStd = BStdSe;
+                    seDev = BDevSe;
+                    stimName = 'B';
+                    switch inclCtrl
+                        case {1,2}
+                            grAvCtrl = BGrAvCtrl;
+                            avCtrl = BAvCtrl;
+                            seCtrl = BSeCtrl;
+                        case 3
+                            grAv50Per = BGrAv50Per;
+                            av50Per = BAv50Per;
+                            se50Per = BSe50Per;
+                            grAvMR = BGrAvMR;
+                            avMR = BAvMR;
+                            seMR = BSeMR;
+                    end
+            end
+            switch plotMtrxI
+                case 0
+                    figure('NumberTitle','off','Name',[combName{c},'_',stimName],'Position',[0,0,1000,600])
+                    title(combNameS(c,s))
+                case 1
+                    switch s % mirror coordinats with s
+                        case 1
+                            nexttile(tileMtrx(combCoord(c-1,2),combCoord(c-1,1)))
+                        case 2
+                            nexttile(tileMtrx(combCoord(c-1,1),combCoord(c-1,2)))
+                    end
+            end
+            switch avType
+                case 1 % grand averages
+                    hold on
+                    switch inclCtrl
+                        case {1,2}
+                            switch plotDifI
+                                case 0 % plot control response only when difference curve is deactivated
+                                    % plot standard error
+                                    seHiCtrl = grAvCtrl(timeCut,filt)+seCtrl(timeCut,filt);
+                                    seLoCtrl = grAvCtrl(timeCut,filt)-seCtrl(timeCut,filt);
+                                    seX = [timetr(timeCut),fliplr(timetr(timeCut))];
+                                    seY = [seHiCtrl',fliplr(seLoCtrl')];
+                                    patch(seX,seY,colCtrl,'FaceAlpha',0.2,'EdgeColor','none')
+                                    % plot data
+                                    plotCtrl = plot(timetr(timeCut),grAvCtrl(timeCut,filt),'color',colCtrl,'Linewidth',2);
+                            end
+                        case 3
+                            switch plotDifI
+                                case 0 % plot control response only when difference curve is deactivated
+                                    % plot standard error
+                                    seHi50Per = grAv50Per(timeCut,filt)+se50Per(timeCut,filt);
+                                    seLo50Per = grAv50Per(timeCut,filt)-se50Per(timeCut,filt);
+                                    seX = [timetr(timeCut),fliplr(timetr(timeCut))];
+                                    seY = [seHi50Per',fliplr(seLo50Per')];
+                                    patch(seX,seY,col50Per,'FaceAlpha',0.2,'EdgeColor','none')
+                                    % plot data
+                                    plot50Per = plot(timetr(timeCut),grAv50Per(timeCut,filt),'color',col50Per,'Linewidth',2);
+                                    % plot standard error
+                                    seHiMR = grAvMR(timeCut,filt)+seMR(timeCut,filt);
+                                    seLoMR = grAvMR(timeCut,filt)-seMR(timeCut,filt);
+                                    seX = [timetr(timeCut),fliplr(timetr(timeCut))];
+                                    seY = [seHiMR',fliplr(seLoMR')];
+                                    patch(seX,seY,colMR,'FaceAlpha',0.2,'EdgeColor','none')
+                                    % plot data
+                                    plotMR = plot(timetr(timeCut),grAvMR(timeCut,filt),'color',colMR,'Linewidth',2);
+                            end
+                    end
+                    % plot standard error
+                    seHiDev = grAvDev(timeCut,filt)+seDev(timeCut,filt);
+                    seLoDev = grAvDev(timeCut,filt)-seDev(timeCut,filt);
+                    seX = [timetr(timeCut),fliplr(timetr(timeCut))];
+                    seY = [seHiDev',fliplr(seLoDev')];
+                    patch(seX,seY,colDev,'FaceAlpha',0.2,'EdgeColor','none')
+                    % plot data
+                    plotDev = plot(timetr(timeCut),grAvDev(timeCut,filt),'color',colDev,'Linewidth',2);
+                    % plot standard error
+                    seHiStd = grAvStd(timeCut,filt)+seStd(timeCut,filt);
+                    seLoStd = grAvStd(timeCut,filt)-seStd(timeCut,filt);
+                    seX = [timetr(timeCut),fliplr(timetr(timeCut))];
+                    seY = [seHiStd',fliplr(seLoStd')];
+                    patch(seX,seY,colStd,'FaceAlpha',0.2,'EdgeColor','none')
+                    % plot data
+                    plotStd = plot(timetr(timeCut),grAvStd(timeCut,filt),'color',colStd,'Linewidth',2);
+                    xticks(timeWin(1):5:timeWin(2))
+                    set(gca,'XTickLabel',timeStart:5:timeEnd)
+                    switch stimMark
+                        case 1 % plot bar
+                            y_dist = y_lim(2)-y_lim(1);
+                            y_min = y_lim(1)+y_dist*0.05;
+                            x_stim = [stimWin(1);stimWin(2);stimWin(2);stimWin(1)];
+                            y_stim = [y_min+y_dist*0.05;y_min+y_dist*0.05;y_min+y_dist*0.07;y_min+y_dist*0.07];
+                            patch('XData',x_stim,'YData',y_stim,'EdgeColor','Black','EdgeAlpha',1,'FaceColor','Black','FaceAlpha',1)
+                    end
+                    switch plotMtrxI
+                        case 0
+                            xlabel('Time [ms]','FontWeight','bold','FontSize',20)
+                            ylabel(ylableName,'FontWeight','bold','FontSize',20)
+                            switch inclCtrl
+                                case 0
+                                    legend([plotDev,plotStd],nameDev,nameStd,'Location','northwest','AutoUpdate','off')
+                                case {1,2}
+                                    switch plotDifI
+                                        case 0
+                                            legend([plotDev,plotStd,plotCtrl],nameDev,nameStd,nameCtrl,'Location','northwest','AutoUpdate','off')
+                                        case 1
+                                            legend([plotDev,plotStd],nameDev,nameStd,'Location','northwest','AutoUpdate','off')
+                                    end
+                                case 3
+                                    switch plotDifI
+                                        case 0
+                                            legend([plotDev,plotStd,plot50Per,plotMR],nameDev,nameStd,name50Per,nameMR,'Location','northwest','AutoUpdate','off')
+                                        case 1
+                                            legend([plotDev,plotStd],nameDev,nameStd,'Location','northwest','AutoUpdate','off')
+                                    end
+                            end
+                        case 1
+                            xlabel(t,'Time [ms]','FontWeight','bold','FontSize',20)
+                            ylabel(t,ylableName,'FontWeight','bold','FontSize',20)
+                    end
+                    set(gca,'FontSize',20,'Linewidth',2,'FontName','Arial')
+                    switch plotWinI
+                        case 1
+                            x_min = timetr(winAll(:,1));
+                            x_max = timetr(winAll(:,2));
+                            x = [x_min;x_max;x_max;x_min];
+                            switch fixedY
+                                case 0
+                                    y_limAuto = ylim;
+                                    y_min = y_limAuto(1);
+                                    y_max = y_limAuto(2);
+                                case 1
+                                    y_min = y_lim(1);
+                                    y_max = y_lim(2);
+                            end
+                            y_temp = [y_min;y_min;y_max;y_max];
+                            y = repmat(y_temp,1,uWin(end));
+                            for w = uWin
+                                switch inclCtrl
+                                    case 0
+                                        statsV = pV(c,s,filt,w);
+                                    case {1,2,3}
+                                        switch statsType
+                                            case 1 % use t-test statistics
+                                                statsV = pV(c,s,filt,w); 
+                                            case 2 % use ANOVA statistics
+                                                statsV = multcompV{c,s,filt,w}.pValue(1); % extract pValue for comparison between dev and std condition
+                                        end
+                                end
+                                if statsV<0.05
+                                    sigWin = patch('XData',x(:,w),'YData',y(:,w),'EdgeColor','black','EdgeAlpha',1,'FaceColor','black','FaceAlpha',.2,'Linewidth',2);
+                                else
+                                    sigWin = patch('XData',x(:,w),'YData',y(:,w),'EdgeColor','black','EdgeAlpha',1,'FaceColor','none','Linewidth',2);
+                                end
+                                uistack(sigWin,'bottom')
+                            end
+                    end
+                    xlim([timeWin(1),timeWin(2)])
+                    if fixedY==1
+                        ylim(y_lim)
+                    end
+                case 2 % individual averages
+                    for f = 1:nFiles
+                        figure ('NumberTitle','off','Name',[combName{c},'_',stimName,'_',convertStringsToChars(recID(f))],'Position',[0,0,1000,500])
+                        hold on
+                        switch inclCtrl
+                            case {1,2}
+                                switch plotDifI
+                                    case 0 % plot control response only when difference curve is deactivated
+                                        plotCtrl = plot(timetr(timeCut),avCtrl(timeCut,f,filt),'color',colCtrl,'Linewidth',2);
+                                end
+                            case 3
+                                switch plotDifI
+                                    case 0 % plot control response only when difference curve is deactivated
+                                        plot50Per = plot(timetr(timeCut),av50Per(timeCut,f,filt),'color',col50Per,'Linewidth',2);
+                                        plotMR = plot(timetr(timeCut),avMR(timeCut,f,filt),'color',colMR,'Linewidth',2);
+                                end
+                        end
+                        plotStd = plot(timetr(timeCut),avStd(timeCut,f,filt),'color',colStd,'Linewidth',2);
+                        plotDev = plot(timetr(timeCut),avDev(timeCut,f,filt),'color',colDev,'Linewidth',2);
+                        xticks(timeWin(1):5:timeWin(2))
+                        set(gca,'XTickLabel',timeStart:5:timeEnd)
+                        switch stimMark
+                            case 1 % plot bar
+                                y_dist = y_lim(2)-y_lim(1);
+                                y_min = y_lim(1)+y_dist*0.05;
+                                x_stim = [stimWin(1);stimWin(2);stimWin(2);stimWin(1)];
+                                y_stim = [y_min+y_dist*0.05;y_min+y_dist*0.05;y_min+y_dist*0.07;y_min+y_dist*0.07];
+                                patch('XData',x_stim,'YData',y_stim,'EdgeColor','Black','EdgeAlpha',1,'FaceColor','Black','FaceAlpha',1)
+                        end
+                        xlabel('Time [ms]','FontWeight','bold','FontSize',20)
+                        ylabel(ylableName,'FontWeight','bold','FontSize',20)
+                        set(gca,'FontSize',20,'Linewidth',2,'FontName','Arial')
+                        xlim([timeWin(1),timeWin(2)])
+                        if fixedY==1
+                            ylim(y_lim)
+                        end
+                    end
+            end
+            if saveI==1
+                switch plotMtrxI
+                    case 0
+                        saveas(gcf,[combName{c},'_',stimName,'_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.jpg'])
+                        saveas(gcf,[combName{c},'_',stimName,'_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.svg'])
+                end
+            end
+        end
+    end
+    if saveI==1
+        switch plotMtrxI
+            case 1
+                saveas(gcf,['Matrix_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.jpg'])
+                saveas(gcf,['Matrix_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.svg'])
+        end
+    end
+end
+%% plotting: data-heatmap
+
+if plotHmI==1
+    % do some prior calculations
+    rgb1 = [(0:0.01:1),ones(1,101)]';
+    rgb2 = [(0:0.01:1),(1:-0.01:0)]';
+    rgb3 = [ones(1,101),(1:-0.01:0)]';
+    cMap = [rgb1,rgb2,rgb3];
+    % create matrix containing the data of interest
+    hmMtrx = NaN(nStim,nStim);
+
+    for c = cStart:nComb
+        for s = 1:2
+            for w = uWin
+                switch hmDataI
+                    case 1
+                        data = round(cohensD_V(c,s,filt,w,1),2); % Cohen's D of comparison between dev and std responses
+                    case 2
+                        data = round(rmsGrAv{c}(s,3,filt,w),2); % grand average RMS values of control response
+                end
+                switch s % mirror coordinats with s
+                    case 1
+                        hmMtrx(tileMtrx(combCoord(c-1,1),combCoord(c-1,2))) = data; % ATTENTION: coordinates change with s in the opposite way as when creating a tiled layout because matrix is filled up from top to bottom while tiled layout is filled from left to right
+                    case 2
+                        hmMtrx(tileMtrx(combCoord(c-1,2),combCoord(c-1,1))) = data;
+                end
+            end
+        end
+    end
+    % create heatmap
+    xvalues = stimOrd;
+    yvalue = stimOrd;
+    figure('NumberTitle','off','Name','Heatmap','Position',[0,0,1450,1300],'Renderer','painters')
+    h = heatmap(xvalues,yvalue,hmMtrx);
+    % do calculations for heatmap
+    cmap = colormap(h);
+    sCmap = size(cmap);
+    h.Title = 'Cohens D';
+    h.XLabel = 'Active Stimulus';
+    h.YLabel = 'Modulatory Stimulus';
+    h.MissingDataColor = [0.7,0.7,0.7];
+    lims = clim;
+    hmMax = max(abs(hmMtrx),[],'all');
+    switch hmDataI
+        case 1
+            hmMin = -hmMax;
+            h.Colormap = cMap;
+        case 2
+            hmMin = 0;
+            h.Colormap = summer;
+    end
+    clim([hmMin,hmMax])
+    set(gca,'FontSize',20,'FontName','Arial')
+    if saveI==1
+        saveas(gcf,['Heatmap_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.jpg'])
+        saveas(gcf,['Heatmap_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.svg'])
+    end
+end
+%% plotting: boxplots
+
+if plotBpI==1
+    switch plotMtrxI
+        case 1
+            figure('NumberTitle','off','Name','Boxpl_All Stims','Position',[0,0,1500,1300],'Renderer','painters')
+            t = tiledlayout(nStim,nStim);
+            title(t,'RMS','FontWeight','bold','FontSize',25)
+    end
+    switch inclCtrl
+        case 0
+            xDev = [1,4,7,10,13];
+            xStd = [2,5,8,11,14];
+            sig_groups = {[1,2],[4,5],[7,8],[10,11],[13,14]};
+        case {1,2}
+            xDev = [1,5,9,13,17];
+            xStd = [2,6,10,14,18];
+            xCtrl = [3,7,11,15,19];
+            sig_groups1 = {[1,2],[5,6],[9,10],[13,14],[17,18]}; % dev vs. std
+            sig_groups2 = {[1,3],[5,7],[9,11],[13,15],[17,19]}; % dev vs. ctrl
+            sig_groups3 = {[2,3],[6,7],[10,11],[14,15],[18,19]}; % std vs. ctrl
+        case 3
+            xDev = [1,6,11,16,21];
+            xStd = [2,7,12,17,22];
+            x50Per = [3,8,13,18,23];
+            xMR = [4,9,14,19,24];
+            sig_groups1 = {[1,2],[6,7],[11,12],[16,17],[21,22]}; % dev vs. std
+            sig_groups2 = {[1,3],[6,8],[11,13],[16,18],[21,23]}; % dev vs. 50Per
+            sig_groups3 = {[1,4],[6,9],[11,14],[16,19],[21,24]}; % dev vs. MR
+            sig_groups4 = {[2,3],[7,8],[12,13],[17,18],[22,23]}; % 50Per vs. std
+            sig_groups5 = {[2,4],[7,9],[12,14],[17,19],[22,24]}; % MR vs. std
+            sig_groups6 = {[3,4],[8,9],[13,14],[18,19],[23,24]}; % MR vs. 50Per
+    end
+
+    for c = cStart:nComb % run once for each stimulus combination
+        for s = 1:2 % once for A, once for B
+            switch s
+                case 1
+                    stimName = 'A';
+                case 2
+                    stimName = 'B';
+            end
+            rmsMin = floor(min(rmsV{c}(:,s,:,filt,:),[],'all')); % identify minimum RMS value for y-limit
+            rmsMax = ceil(max(rmsV{c}(:,s,:,filt,:),[],'all')); % identify maximum RMS value for y-limit
+            switch plotMtrxI
+                case 0
+                    figure('NumberTitle','off','Name',['Boxpl_',combName{c},'_',stimName],'Position',[0,0,boxpFigWi,400])
+                    title('RMS')
+                case 1
+                    switch s % mirror coordinats with s
+                        case 1
+                            nexttile(tileMtrx(combCoord(c-1,2),combCoord(c-1,1)))
+                        case 2
+                            nexttile(tileMtrx(combCoord(c-1,1),combCoord(c-1,2)))
+                    end
+            end
+            ylim([rmsMin-0.5,rmsMax+1])
+            for w = uWin
+                hold on
+                boxDev = boxchart(xDev(w)*ones(size(rmsV{c}(:,s,1,filt,w))),rmsV{c}(:,s,1,filt,w),'BoxFaceColor',colDev,'MarkerColor',colDev,'Linewidth',2);
+                boxStd = boxchart(xStd(w)*ones(size(rmsV{c}(:,s,2,filt,w))),rmsV{c}(:,s,2,filt,w),'BoxFaceColor',colStd,'MarkerColor',colStd,'Linewidth',2);
+                % switch inclCtrl
+                %     case 0
+                %         H = sigstar(sig_groups{w},pV(c,s,filt,w),cohensD_V(c,s,filt,w,1));
+                %     case {1,2}
+                %         boxCtrl = boxchart(xCtrl(w)*ones(size(rmsV{c}(:,s,3,filt,w))),rmsV{c}(:,s,3,filt,w),'BoxFaceColor',colCtrl,'MarkerColor',colCtrl,'Linewidth',2);
+                %         switch statsType
+                %             case 1
+                %                 H = sigstar(sig_groups1{w},pV(c,s,filt,w),cohensD_V(c,s,filt,w,1));
+                %             case 2
+                %                 statsV = multcompV{c,s,filt,w}.pValue([1,2,4]);
+                %                 eS = zeros(3,1);
+                %                 for cmp = 1:3
+                %                     eS(cmp) = cohensD_V(c,s,filt,w,cmp);
+                %                 end
+                %                 H = sigstar({sig_groups1{w},sig_groups2{w},sig_groups3{w}},statsV,eS);
+                %         end
+                %     case 3
+                %         box50Per = boxchart(x50Per(w)*ones(size(rmsV{c}(:,s,3,filt,w))),rmsV{c}(:,s,3,filt,w),'BoxFaceColor',col50Per,'MarkerColor',col50Per,'Linewidth',2);
+                %         boxMR = boxchart(xMR(w)*ones(size(rmsV{c}(:,s,4,filt,w))),rmsV{c}(:,s,4,filt,w),'BoxFaceColor',colMR,'MarkerColor',colMR,'Linewidth',2);
+                %         switch statsType
+                %             case 1
+                %                 H = sigstar(sig_groups1{w},pV(c,s,filt,w),cohensD_V(c,s,filt,w,1));
+                %             case 2
+                %                 statsV = multcompV{c,s,filt,w}.pValue([1,2,3,5,6,9]);
+                %                 eS = zeros(6,1);
+                %                 for cmp = 1:6
+                %                     eS(cmp) = cohensD_V(c,s,filt,w,cmp);
+                %                 end
+                %                 H = sigstar({sig_groups1{w},sig_groups2{w},sig_groups3{w},sig_groups4{w},sig_groups5{w},sig_groups6{w}},statsV,eS);
+                %         end
+                % end
+                set(H(:,2),'FontSize',25,'FontName','Arial')
+            end
+            switch inclCtrl
+                case 0
+                    xlim([xDev(uWin(1))-1,xStd(uWin(end))+1])
+                    xticks((xDev(uWin)+xStd(uWin))/2)
+%                     legend([boxDev,boxStd],nameDev,nameStd,'Location','northwest','AutoUpdate','off')
+                case {1,2}
+                    xlim([xDev(uWin(1))-1,xCtrl(uWin(end))+1])
+                    xticks(xStd(uWin))
+%                     legend([boxDev,boxStd,boxCtrl],nameDev,nameStd,nameCtrl,'Location','northwest','AutoUpdate','off')
+                case 3
+                    xlim([xDev(uWin(1))-1,xMR(uWin(end))+1])
+                    xticks((xStd(uWin)+x50Per(uWin))/2)
+%                     legend([boxDev,boxStd,box50Per,boxMR],nameDev,nameStd,name50Per,nameMR,'Location','northwest','AutoUpdate','off')
+            end
+            switch plotMtrxI
+                case 0
+                    ylabel(ylableName,'FontWeight','bold','FontSize',20)
+                    xticklabels(winTxt(uWin))
+                case 1
+                    ylabel(t,ylableName,'FontWeight','bold','FontSize',20)
+                    xticks([])
+                    ylim([0,15])
+            end
+            set(gca,'FontSize',25,'Linewidth',2,'FontName','Arial')
+            if saveI==1
+                switch plotMtrxI
+                    case 0
+                        saveas(gcf,['Boxp_',combName{c},'_',stimName,'_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.jpg'])
+                        saveas(gcf,['Boxp_',combName{c},'_',stimName,'_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.svg'])
+                end
+            end
+        end
+    end
+    if saveI==1
+        switch plotMtrxI
+            case 1
+                saveas(gcf,['Boxpl_Matrix_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.jpg'])
+                saveas(gcf,['Boxpl_Matrix_',filtName,'_cutTime',num2str(cutTime),'_ctrl',num2str(inclCtrl),'.svg'])
+        end
+    end
+end