123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338 |
- %Figure 3D: boxplots of change in node strenght relative to baseline
- %get data
- cd 'D:\Labor\Analyses\Toolbox\fMRI'
- addpath('FSLNets')
- addpath('L1precision')
- addpath('pwling')
- addpath('../Tools/BCT')
- folder_idx = pwd; % get folder ID to store data
- files_idx = dir([folder_idx,'/*mat']); % list data2load
- baseFileName = files_idx.name;
- fullFileName = fullfile(folder_idx, baseFileName);
- dirInfos = dir(fullFileName);
- clear ('folder_idx','files_idx','baseFileName','fullFileName');
- cd 'D:\Labor\Projects\Stroke_tDCS_rsfmri\B_rsfmri_stroke_tDCS_mice\data';
- targetFolder=pwd;
- clear 'data';
- load data.mat;
- load labels.mat;
- group_name = {'Control','PT','PT_tDCS'};
- %% get individual node values
- tick = 1;
- for i = 1:size(data,2)
- group = data(i).group;
- for k = 1:size(data(i).subj,2)
- sbj = data(i).subj(k).ID;
- for t = 1:size(data(i).subj(k).tp,2)
- mat = data(i).subj(k).tp(t).mNet;
- if ~isempty(mat)
- for node = 1:length(labels.nr14)
-
- node_strength = sum(mat(labels.nr14(node),labels.nr14),'omitnan');
- master{tick,node+3} = node_strength;
-
- end
- master{tick,1} = group;
- master{tick,2} = sbj;
- master{tick,3} = num2str(t);
- tick = tick+1;
- end
- end
- end
- end
- header = {'group','subject','timepoint'};
- for kx = 1:length(labels.ROI14_short)
- header(kx+3)=labels.ROI14_short(kx);
- end
- master_tab = cell2table(master,'VariableNames',header);
- master_tab.group = categorical(master_tab.group);
- master_tab.subject = categorical(master_tab.subject);
- master_tab.timepoint = categorical(master_tab.timepoint);
- %%
- tack = 1;
- for ix = 1:length(group_name)
- subdata= master_tab(master_tab.group == group_name{ix},:);
- IDlist = unique(subdata.subject);
- % get group mean baseline values
- for ox = 1:length(labels.nr14)
- baseValues(ox) = mean(cell2mat(table2cell(subdata(subdata.timepoint=='1',ox+3))));
- end
- % loop over single subjects
- for mx = 1:length(IDlist)
- %determine if subject data at timepoint is available
- if sum(subdata.subject==IDlist(mx) & subdata.timepoint=='3')~=0
- IDdata = cell2mat(table2cell(subdata(subdata.subject==IDlist(mx) & subdata.timepoint=='3',4:end)));
-
- % determine if subject has individual baseline
- if sum(subdata.subject==IDlist(mx) & subdata.timepoint=='1')~= 0
- baseData = cell2mat(table2cell(subdata(subdata.subject==IDlist(mx) & subdata.timepoint=='1',4:end)));
- deltaValues = IDdata - baseData;
- else
- deltaValues = IDdata - baseValues;
- end
- Delta{tack,1} = group_name{ix};
- Delta{tack,2} = string(IDlist(mx));
- for lx = 1:length(labels.ROI14_short)
- Delta{tack,lx+3} = deltaValues(lx);
- end
- tack = tack+1;
- end
- end
- end
- %% get control shift
- for ox = 1:length(labels.nr14)
-
- baseValues(1,ox) = mean(cell2mat(table2cell(...
- master_tab(master_tab.group == 'Control' & ...
- master_tab.timepoint =='1',ox+3)...
- )));
- baseValues(2,ox) = mean(cell2mat(table2cell(...
- master_tab(master_tab.group == 'Control' & ...
- master_tab.timepoint =='3',ox+3)...
- )));
- end
- Delta_base = baseValues(2,:)-baseValues(1,:);
- % Delta(:,4:end) = cell2mat(Delta(:,4:end)) - Delta_base;
- %%
- % get concanated list of data
- for kx = 1:length(labels.ROI14_short)
- subseg = Delta(:,[1:3 kx+3]);
- DeltaList([1+(kx-1)*size(subseg,1):(kx*size(subseg,1))],:)= subseg;
- DeltaList([1+(kx-1)*size(subseg,1):(kx*size(subseg,1))],3)= labels.ROI14_short(kx);
- end
- % transform cell to Table
- DeltaHead = {'group','subject','ROI','value'};
- Delta_Tab = cell2table(DeltaList,'VariableNames',DeltaHead);
- Delta_Tab.ROI = categorical(Delta_Tab.ROI);
- bpColors (1,:) = hex2rgb('777777');
- bpColors (2,:) = hex2rgb('D95319');
- bpColors (3,:) = hex2rgb('0072BD');
- %%
- PT_Delta = Delta;
- GroupID = categorical(PT_Delta(:,1));
- rawdata = cell2mat(PT_Delta(GroupID=='PT_tDCS',3:end));
- for r = 1:size(rawdata,2)
- pd = fitdist(rawdata(:,r),'Normal');
- ci = paramci(pd);
- statPT(4,r) = pd.mean;
- statPT(5:6,r) = ci(:,1);
- end
- statPT ([1 4],:) = statPT([1 4],:)- Delta_base;
- tabNodePT = array2table(statPT,'VariableNames',labels.ROI14_short)
- %% Permutaion Based Statistics (PBS) of node strength
- eval = squeeze(mean(cell2mat(PT_Delta(GroupID=='PT',3:end)))) - ...
- squeeze(mean(cell2mat(PT_Delta(GroupID=='PT_tDCS',3:end)))) ;
- % number of permutations
- n_permutes = 1000;
- conc_data= cell2mat(PT_Delta(:,3:end));
- nx=2;%*size(conc_data,2); % Anzahl an tests; Bonferroni Korrektur?
- zval_5 = abs(norminv(0.05/nx)); % for p=0.05
- zval_1 = abs(norminv(0.01/nx)); % for p=0.01
- zval_01 = abs(norminv(0.001/nx)); % for p=0.01
- zval=abs(norminv(0.025));
- for i=50000:-1:1
- p=i/100000;
- pval(i)=abs(norminv(p));
- end
- for permi=1:n_permutes
- randorder = randperm(size(conc_data,1));
- temp = conc_data(randorder,:);
- perm_map(permi,:) = squeeze(mean(temp(1:floor(size(conc_data,1)/2),:))) - ...
- squeeze(mean(temp(ceil(size(conc_data,1)/2):end,:)));
- end
- mean_h0 = squeeze(mean(perm_map,1));
- std_h0 = std(perm_map,1);
- zstat = abs(eval-mean_h0)./(std_h0);
- idx=zstat>zval_5;
- imx=zstat>zval;
- mat=eval.*idx;
- for i=1:size(zstat,1)
- for k=1:size(zstat,2)
- [c idx] = min(abs(pval-zstat(i,k)));
- pstat(i,k)=idx/100000;
- end
- end
- FDR=5;
- fr=FDR/100;
- [hidx, crit_p, adj_ci_cvrg, adj_p]=fdr_bh(pstat,fr,'dep');
- hidx
- %%
- [B,I] = sort(eval,'descend');
- ROI_order = labels.ROI14_short(I);
- Delta_Tab.ROI = categorical(Delta_Tab.ROI,ROI_order);
- figure
- b = boxchart(Delta_Tab.ROI,Delta_Tab.value,'GroupByColor',Delta_Tab.group,'BoxWidth',0.75);
- for bx=1:3
- b(bx).BoxFaceColor = bpColors(bx,:);
- b(bx).MarkerColor = bpColors(bx,:);
- end
- ylabel('Delta Node Strength (Day 14)');
- legend({'Control','PT sham','PT tDCS'})
- ax = gca;
- ax.FontName = 'Times';
- ax.FontSize = 22;
- %%
- %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- %Figure 3E-F: Group Mean Boxplots of sensorimotor network
- %%
- master_tab = DataPTrsfmri;
- group_name = {'control','PT sham','PT tDCS'};
- %% get control shift
- for ox = 1:length(group_name)
-
- baseValues(1,:) = squeeze(mean(cell2mat(table2cell(...
- master_tab(master_tab.Group == group_name(ox) & ...
- master_tab.Time =='baseline',3:end)...
- )),1,'omitnan'));
- baseValues(2,:) = squeeze(mean(cell2mat(table2cell(...
- master_tab(master_tab.Group == group_name(ox) & ...
- master_tab.Time =='Day 28',3:end)...
- )),1,'omitnan'));
- Delta_base(ox,:) = baseValues(2,:)-baseValues(1,:);
- end
- % Delta(:,4:end) = cell2mat(Delta(:,4:end)) - Delta_base;
- %%
- rawdata = cell2mat(table2cell(master_tab(master_tab.Group=='control'...
- & master_tab.Time=="Day 28",3:end)));
- for r = 1:size(rawdata,2)
- pd = fitdist(rawdata(:,r),'Normal');
- ci = paramci(pd);
- stat(1,r) = pd.mean;
- stat([2 3],r) = ci(:,1);
- end
- % statPT ([1 4],:) = statPT([1 4],:)- Delta_base;
- %
- % tabNodePT = array2table(statPT,'VariableNames',labels.ROI14_short)
- %%
- %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- %Figure 3A-D: Group Mean Boxplots of sensorimotor network
- %get data
- cd 'D:\Labor\Analyses\Toolbox\fMRI'
- addpath('FSLNets')
- addpath('L1precision')
- addpath('pwling')
- addpath('../Tools/BCT')
- folder_idx = pwd; % get folder ID to store data
- files_idx = dir([folder_idx,'/*mat']); % list data2load
- baseFileName = files_idx.name;
- fullFileName = fullfile(folder_idx, baseFileName);
- dirInfos = dir(fullFileName);
- clear ('folder_idx','files_idx','baseFileName','fullFileName');
- cd 'D:\Labor\Projects\Stroke_tDCS_rsfmri\B_rsfmri_stroke_tDCS_mice\data';
- targetFolder=pwd;
- clear 'data';
- load data.mat;
- load labels.mat;
- group_name = {'Control','PT','PT_tDCS'};
- %%
- for ri = 1:size(data,2)
- for li = 1:size(data(ri).subj,2)
- for ti = 1:size(data(ri).subj(li).tp,2)
-
- mnx=data(ri).subj(li).tp(ti).mNet;
- if ~isempty(mnx)
- mNet(:,:,(ri-1)*4+ti,li)=mnx;
- end
-
- end
- % rnx=data(ri).subj(li).mNet;
- % ROI(:,:,li,ri)=rnx;
- end
- end
- mNet(mNet==0)=nan;
- clear ('mnx','rnx','li','ri');
- mNet_mean = squeeze(mean(mNet,4,'omitnan'));
- for i = 1:3
- for k=1:4
- DeltaBaseline (:,:,k+(i-1)*4) = mNet_mean(:,:,k+(i-1)*4)-mNet_mean(:,:,1+(i-1)*4);
- end
- end
- ROI_Delta = DeltaBaseline(labels.nr14,labels.nr14,:);
- %%
- for i = 1:size(ROI_Delta,3)
- plotmat_values_rsfMRI(ROI_Delta(:,:,i),targetFolder,dirInfos,[-.75 .75],0,labels.ROI14_short);
- hold on
- axis square;
- plot(get(gca,'xlim'),[mean(get(gca,'ylim')), mean(get(gca,'ylim'))],'k-','linewidth',2);
- plot([mean(get(gca,'xlim')), mean(get(gca,'xlim'))], get(gca,'ylim'),'k-','linewidth',2);
- title(num2str(i))
- % title(strcat(name))%,'; mean node strength = ', num2str(norm_ROI_delta(ix*4+k))));
- ax=gca;
- ax.FontSize=20;
- ax.FontName='Times';
- end
|