Skip to content
Snippets Groups Projects
Commit 54a195ad authored by Tamas Borbath's avatar Tamas Borbath
Browse files

Cleaned up a bit and adapted plot_all_test_data_profit code to new... 

Cleaned up a bit and adapted plot_all_test_data_profit code to new Bland-Altman plots. However the code is still not clean
parent 92689f4a
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
createTestData/plot_all_test_data_profit.m
+ 183
− 58
View file @ 54a195ad
......@@ -27,12 +27,8 @@ exportFolders = {'0 Water References',...
pathNameExport = 'ProFit test data';
pathBaseExportFiles = pathToDataFolder(pathNameExport, exportFolders);
% pathBaseExportFiles = pathBaseExportFiles(1:end-1);
% pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 trunc\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 sinebell\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 matched\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - R1 0.25scale m_15 trunc\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - R2 0.25scale m_15 trunc\'];
pathBaseExportFiles = pathBaseExportFiles(1:end-1);
pathBaseExportFiles = [pathBaseExportFiles ' - R3\'];
dataExportPathBaseRef = strcat(pathBaseExportFiles, exportFolders{2}, '\');
dataExportPathBase = strcat(pathBaseExportFiles, exportFolders{reconOption+1}, '\');
......@@ -105,7 +101,7 @@ if fitLCModel
else
SoftwareName = 'ProFit';
end
[metConcAllSubjLCModel, metConcAllSubjRefLCModel, metConcAllSubjNormLCModel] = ...
[metConcAllSubjLCModel, metConcAllSubjRefLCModel, metConcAllSubjNormLCModel, MC_H2O_ConcAllSubjLCModel, MC_H2O_ConcAllSubjRefLCModel] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
......@@ -117,7 +113,7 @@ if fitLCModel
else
SoftwareName = 'ProFit';
end
[metConcAllSubjProFit, metConcAllSubjRefProFit, metConcAllSubjNormProFit] = ...
[metConcAllSubjProFit, metConcAllSubjRefProFit, metConcAllSubjNormProFit, MC_H2O_ConcAllSubjProFit, MC_H2O_ConcAllSubjRefProFit] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
......@@ -136,15 +132,17 @@ metConcAllSubjLCModelScaled = metConcAllSubjLCModel;
metConcAllSubjRefLCModelScaled = metConcAllSubjRefLCModel;
% BA plot paramters
tit = 'Fitting Repeatability '; % figure title
subsets = {'32 ave.', '64 ave.'};
% subsets = {'32 ave.1','32 ave.2','32 ave.3', '64 ave.1', '64 ave.2', '64 ave.3'};
tit = 'Repeatability '; % figure title
subsets = {'32 ave.1','32 ave.2','32 ave.3', '64 ave.1', '64 ave.2', '64 ave.3'};
subsets_refit = {'32 ave.', '64 ave.'};
softwareNames = {'ProFit', 'LCModel'};
gnames = {subsets, softwareNames}; % names of groups in data {dimension 1 and 2}
label = {'Ref. Fits','Subset Fits','ratio'}; % Names of data sets
label_refit = {'c^{fits1}_{i,k}','c^{fits2}_{i,k}','arb.u.'}; % Names of data sets
label_SW_comp = {'LCModel','ProFit','arb.u.'}; % Names of data sets
corrinfo = {'n','SSE','r2','eq'}; % stats to display of correlation scatter plot
BAinfo = {'RPC','ks'}; % stats to display on Bland-ALtman plot
BAinfo_refit = {''}; % stats to display on Bland-ALtman plot
limits = 'auto'; % how to set the axes limits
if 1 % colors for the data sets may be set as:
colors = 'br'; % character codes
......@@ -152,61 +150,147 @@ else
colors = [0 0 1;... % or RGB triplets
1 0 0];
end
paramsBA.data1TreatmentMode = 'Compare';
paramsBA.diffValueMode = 'percent';
numOfSimulations = length(namingSuffixesAve)*length(namingSuffixesCoil);
rpc_table = cell(numberOfMet+3,3);
%there seems to be a fit error in the fitting of even the MC water data.
% BlandAltman(MC_H2O_ConcAllSubjLCModel, MC_H2O_ConcAllSubjProFit*1e-6, label_SW_comp,[tit 'Water'],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'showFitCI',' on',...
% 'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
%
% BlandAltman(MC_H2O_ConcAllSubjRefLCModel, MC_H2O_ConcAllSubjRefProFit*1e-6, label_SW_comp,[tit 'Water'],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'showFitCI',' on',...
% 'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
rpc_table = cell(numberOfMet+5,5);
rpc_table{1,2} = softwareNames{1};
rpc_table{1,3} = softwareNames{2};
rpc_table(2:end-2,1) = metabolitesToDisplay;
rpc_table{1,4} = softwareNames{2};
rpc_table{2,2} = subsets_refit{1};
rpc_table{2,3} = subsets_refit{2};
rpc_table{2,4} = subsets_refit{1};
rpc_table{2,5} = subsets_refit{2};
rpc_table(3:end-3,1) = metabolitesToDisplay;
rpc_table{end-1,1} = 'Mean Main';
rpc_table{end,1} = 'Mean';
%title('Reproducibility');
subplotsX = 4;
subplotsY = 4;
totalSubplots = subplotsX * subplotsY;
rpcSum = zeros(1,numberOfMet);
diffMedian = zeros(1,numberOfMet);
for indexMet = 1: numberOfMet
if indexMet == metaboliteRef
continue;
end
if mod(indexMet,totalSubplots/2) == 1
figure;
end
runningIndexPlot = (indexMet-1)+1+floor((indexMet-1)/subplotsY)*subplotsY-floor((indexMet-1)/totalSubplots*2)*totalSubplots;
%ProFit Concentrations
metConcAllSubjProFitPerCr_ = metConcAllSubjProFitScaled(:, :, indexMet);
metConcAllSubjProFitPerCr_ = reshape(metConcAllSubjProFitPerCr_,11*3,2);%group same SNR
metConcAllSubjRefProFitPerCr_ = metConcAllSubjRefProFitScaled(:, indexMet);
metConcAllSubjRefProFitPerCr_ = repmat(metConcAllSubjRefProFitPerCr_,1,6);
metConcAllSubjRefProFitPerCr_ = reshape(metConcAllSubjRefProFitPerCr_,11*3,2);%group same SNR
metConcAllSubjProFitScaled_ = metConcAllSubjProFitScaled(:, :, indexMet);
metConcAllSubjProFitScaled_ = reshape(metConcAllSubjProFitScaled_,11*3,2);%group same SNR
metConcAllSubjRefProFitScaled_ = metConcAllSubjRefProFitScaled(:, indexMet);
metConcAllSubjRefProFitScaled_ = repmat(metConcAllSubjRefProFitScaled_,1,6);
metConcAllSubjRefProFitScaled_ = reshape(metConcAllSubjRefProFitScaled_,11*3,2);%group same SNR
%LCModel Concentrations
metConcAllSubjLCModelPerCr_ = metConcAllSubjLCModelScaled(:, :, indexMet);
metConcAllSubjLCModelPerCr_ = reshape(metConcAllSubjLCModelPerCr_,11*3,2);%group same SNR
metConcAllSubjRefLCModelPerCr_ = metConcAllSubjRefLCModelScaled(:, indexMet);
metConcAllSubjRefLCModelPerCr_ = repmat(metConcAllSubjRefLCModelPerCr_,1,6);
metConcAllSubjRefLCModelPerCr_ = reshape(metConcAllSubjRefLCModelPerCr_,11*3,2);%group same SNR
metConcAllSubjLCModelScaled_ = metConcAllSubjLCModelScaled(:, :, indexMet);
metConcAllSubjLCModelScaled_ = reshape(metConcAllSubjLCModelScaled_,11*3,2);%group same SNR
metConcAllSubjRefLCModelScaled_ = metConcAllSubjRefLCModelScaled(:, indexMet);
metConcAllSubjRefLCModelScaled_ = repmat(metConcAllSubjRefLCModelScaled_,1,6);
metConcAllSubjRefLCModelScaled_ = reshape(metConcAllSubjRefLCModelScaled_,11*3,2);%group same SNR
%ProFit&LCModel Concentrations
metConcAllSubjPerCr_ = cat(3, metConcAllSubjProFitPerCr_, metConcAllSubjLCModelPerCr_);
metConcAllSubjRefPerCr_ = cat(3, metConcAllSubjRefProFitPerCr_, metConcAllSubjRefLCModelPerCr_);
%ref Fits vs subfits
% BlandAltman(metConcAllSubjRefPerCr_, metConcAllSubjPerCr_, label,[tit metabolitesLCModel{indexMet}],gnames,...
metConcAllSubjScaled_ = cat(3, metConcAllSubjProFitScaled_, metConcAllSubjLCModelScaled_);
metConcAllSubjRefScaled_ = cat(3, metConcAllSubjRefProFitScaled_, metConcAllSubjRefLCModelScaled_);
%% ref Fits vs subfits: Code not used anymore in this style
% BlandAltman(metConcAllSubjRefScaled_, metConcAllSubjScaled_, label,[tit metabolitesLCModel{indexMet}],gnames,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',colors, 'showFitCI',' on',...
% 'baStatsMode','non-parametric','forceZeroIntercept','on');
%LCModel vs ProFit
BlandAltman(metConcAllSubjLCModelPerCr_, metConcAllSubjProFitPerCr_, label_SW_comp,[tit metabolitesLCModel{indexMet}],subsets,...
'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',colors, 'showFitCI',' on',...
'baStatsMode','non-parametric','forceZeroIntercept','on', 'diffValueMode', 'percent');
%ProFit ref Fits vs subfits
% rpcProFit = BlandAltman(metConcAllSubjRefProFitPerCr_, metConcAllSubjProFitPerCr_, label,[tit softwareNames{1} ' ' metabolitesLCModel{indexMet}],subsets,...
% % Reference fits vs all per Software
% %ProFit ref Fits vs subfits
% rpcProFit = BlandAltman(metConcAllSubjRefProFitScaled_, metConcAllSubjProFitScaled_, label,[tit softwareNames{1} ' ' metabolitesLCModel{indexMet}],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',[0.75 0.5 0; 0.5 0.2 0], 'showFitCI',' on',...
% 'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
% %LCModel ref Fits vs subfits
% rpcLCModel = BlandAltman(metConcAllSubjRefLCModelPerCr_, metConcAllSubjLCModelPerCr_, label,[tit softwareNames{2} ' ' metabolitesLCModel{indexMet}],subsets,...
% rpcLCModel = BlandAltman(metConcAllSubjRefLCModelScaled_, metConcAllSubjLCModelScaled_, label,[tit softwareNames{2} ' ' metabolitesLCModel{indexMet}],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',[0, 0.5, 0.75; 0, 0.25 0.5], 'showFitCI',' on',...
% 'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
% rpc_table{indexMet+1,2} = rpcProFit;
% rpc_table{indexMet+1,3} = rpcLCModel;
% %rpc_table{indexMet+1,2} = rpcProFit;
% %rpc_table{indexMet+1,3} = rpcLCModel;
%% LCModel vs ProFit
% BlandAltman(metConcAllSubjLCModelScaled_, metConcAllSubjProFitScaled_, label_SW_comp,[tit metabolitesLCModel{indexMet}],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',colors, 'showFitCI',' on',...
% 'baStatsMode','non-parametric','forceZeroIntercept','on', 'diffValueMode', 'percent');
% [rpc, ~, stat] = BlandAltman(metConcAllSubjRefLCModelScaled_, metConcAllSubjRefProFitScaled_, label_SW_comp,[tit metabolitesLCModel{indexMet}],subsets,...
% 'corrInfo',corrinfo,'baInfo',BAinfo,'axesLimits',limits,'colors',colors, 'showFitCI',' on',...
% 'baStatsMode','non-parametric','forceZeroIntercept','on', 'diffValueMode', 'percent');
% rpcs(indexMet) = rpc;
% diffMedian(indexMet) = stat.differenceMedian;
%% Subfits vs subfits per Software
%ProFit ref Fits vs subfits
rpc_ProFit_32_ave = calculate_rpc_percent(metConcAllSubjProFitScaled_(1:11,1), metConcAllSubjProFitScaled_(12:22,1));
rpc_ProFit_64_ave = calculate_rpc_percent(metConcAllSubjProFitScaled_(1:11,2), metConcAllSubjProFitScaled_(12:22,2));
legendText = {['{\color{blue}\bf32-ave: RPC=' mynum2str(rpc_ProFit_32_ave,2) '%}'];['{\color{red}\bf64-ave: RPC=' mynum2str(rpc_ProFit_64_ave,2) '%}']};
sProFit = subplot(subplotsX,subplotsY,runningIndexPlot);
[rpcProFit, fig, stats] = BlandAltmanOnly(sProFit, metConcAllSubjProFitScaled_(1:11,:), metConcAllSubjProFitScaled_(12:22,:),...
label_refit,[metabolitesToDisplay{indexMet}],legendText,...
'corrInfo',corrinfo,'baInfo',BAinfo_refit,'axesLimits',limits,'colors','br', 'showFitCI',' on',...
'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
%LCModel ref Fits vs subfits
rpc_LCModel_32_ave = calculate_rpc_percent(metConcAllSubjLCModelScaled_(1:11,1), metConcAllSubjLCModelScaled_(12:22,1));
rpc_LCModel_64_ave = calculate_rpc_percent(metConcAllSubjLCModelScaled_(1:11,2), metConcAllSubjLCModelScaled_(12:22,2));
legendText = {['{\color{blue}\bf32-ave: RPC=' mynum2str(rpc_LCModel_32_ave,2) '%}'];['{\color{red}\bf64-ave: RPC=' mynum2str(rpc_LCModel_64_ave,2) '%}']};
sLCModel = subplot(subplotsX,subplotsY,runningIndexPlot+subplotsY);
[rpcLCModel, fig, stats] = BlandAltmanOnly(sLCModel, metConcAllSubjRefLCModelScaled_(1:11,:), metConcAllSubjLCModelScaled_(12:22,:), ...
label_refit,[metabolitesToDisplay{indexMet}], legendText,...
'corrInfo',corrinfo,'baInfo',BAinfo_refit,'axesLimits',limits,'colors','br', 'showFitCI',' on',...
'baStatsMode','Gaussian','forceZeroIntercept','on', 'diffValueMode', 'percent');
if mod(indexMet,subplotsY)==1
pos1 = sProFit.Position;
annotation('textarrow', [pos1(1), 0], [pos1(2)+pos1(4)/2, 0],'String',softwareNames{1}, ...
'HeadStyle','none','LineStyle', 'none', 'TextRotation',90,'units','normalized', 'FontSize', 16, ...
'FontWeight','bold','HorizontalAlignment','center', 'VerticalAlignment','middle');
pos2 = sLCModel.Position;
annotation('textarrow', [pos1(1), 0], [pos2(2)+pos2(4)/2, 0],'String',softwareNames{2}, ...
'HeadStyle','none','LineStyle', 'none', 'TextRotation',90,'units','normalized', 'FontSize', 16, ...
'FontWeight','bold','HorizontalAlignment','center', 'VerticalAlignment','middle');
end
rpc_table{indexMet+2,2} = rpc_ProFit_32_ave;
rpc_table{indexMet+2,3} = rpc_ProFit_64_ave;
rpc_table{indexMet+2,4} = rpc_LCModel_32_ave;
rpc_table{indexMet+2,5} = rpc_LCModel_64_ave;
end
rpc_meanProFit = mean(cell2mat(rpc_table(2:end-2,2)));
rpc_meanLCModel = mean(cell2mat(rpc_table(2:end-2,3)));
rpc_table{end,2} = rpc_meanProFit;
rpc_table{end,3} = rpc_meanLCModel;
rpc_values = cell2mat(rpc_table(3:end-2,2:5));
rpc_meanProFit_32ave_main = mean(rpc_values(indecesMainMet,1));
rpc_meanProFit_64ave_main = mean(rpc_values(indecesMainMet,2));
rpc_meanLCModel_32ave_main = mean(rpc_values(indecesMainMet,3));
rpc_meanLCModel_64ave_main = mean(rpc_values(indecesMainMet,4));
rpc_table{end-1,2} = rpc_meanProFit_32ave_main;
rpc_table{end-1,3} = rpc_meanProFit_64ave_main;
rpc_table{end-1,4} = rpc_meanLCModel_32ave_main;
rpc_table{end-1,5} = rpc_meanLCModel_64ave_main;
rpc_meanProFit_32ave = mean(rpc_values(:,1));
rpc_meanProFit_64ave = mean(rpc_values(:,2));
rpc_meanLCModel_32ave = mean(rpc_values(:,3));
rpc_meanLCModel_64ave = mean(rpc_values(:,4));
rpc_table{end,2} = rpc_meanProFit_32ave;
rpc_table{end,3} = rpc_meanProFit_64ave;
rpc_table{end,4} = rpc_meanLCModel_32ave;
rpc_table{end,5} = rpc_meanLCModel_64ave;
xlswrite([pathBaseExportFiles, 'InVivoResults_RPC.xlsx'], rpc_table)
%% code will break here. This used to be the old style plotting/results without Bland-Altman plots.
meanMetDeviation = mean(abs(metaboliteConcentrationsAllSubjNorm), 'omitnan');
stdMetDeviation = std(abs(metaboliteConcentrationsAllSubjNorm), 'omitnan');
......@@ -265,7 +349,8 @@ end
% metaboliteConcentrations(1) = metaboliteConcentrations(1) / 1e6; % adjust the MMB to be within range
end
function [metaboliteConcentrationsAllSubj, metaboliteConcentrationsAllSubjRef, metaboliteConcentrationsAllSubjNorm] = ...
function [metaboliteConcentrationsAllSubj, metaboliteConcentrationsAllSubjRef, metaboliteConcentrationsAllSubjNorm,...
MC_H2O_ConcentrationsAllSubj, MC_H2O_ConcentrationsAllSubjRef] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
......@@ -308,25 +393,34 @@ meanFQN = 0;
meanFQN2 = 0;
numOfSimulations = length(namingSuffixesAve)*length(namingSuffixesCoil);
metaboliteConcentrationsAllSubjNorm = zeros(numberOfSubjects * numOfSimulations, numberOfMet);
MC_H2O_ConcentrationsAllSubj = zeros(numberOfSubjects, numOfSimulations);
metaboliteConcentrationsAllSubj = zeros(numberOfSubjects, numOfSimulations, numberOfMet);
metaboliteConcentrationsAllSubjRef = zeros(numberOfSubjects, numberOfMet);
MC_H2O_ConcentrationsAllSubjRef = zeros(numberOfSubjects, 1);
%% do the actual plotting
for indexSubj = 1:numberOfSubjects
if fitLCModel
load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations');
load([dataExportPathRef{indexSubj}, 'concentrations', truncSuffix '_MC_H2O.mat'], 'concentrations_MC_H2O');
concentrationsRef = concentrations;
concentrations_MC_H2ORef = concentrations_MC_H2O;
load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations');
load([dataExportPath{indexSubj}, 'concentrations', truncSuffix '_MC_H2O.mat'], 'concentrations_MC_H2O');
else
load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2');
load([dataExportPathRef{indexSubj}, 'concentrations', truncSuffix '_H2O_profit.mat'], 'concentrations_H2O');
load([dataExportPathRef{indexSubj}, 'concentrations', truncSuffix '_MC_H2O_profit.mat'], 'concentrations_MC_H2O');
concentrationsRef = concentrations;
concentrations_H2ORef = cell2mat(concentrations_H2O);
concentrations_MC_H2ORef = cell2mat(concentrations_MC_H2O);
activeMetabolitesRef = activeMetabolites;
FQNsRef = FQNs;
FQNsRef2 = FQNs2;
load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2');
load([dataExportPath{indexSubj}, 'concentrations', truncSuffix '_H2O_profit.mat'], 'concentrations_H2O');
load([dataExportPath{indexSubj}, 'concentrations', truncSuffix '_MC_H2O_profit.mat'], 'concentrations_MC_H2O');
concentrations_H2O = cell2mat(concentrations_H2O);
concentrations_MC_H2O = cell2mat(concentrations_MC_H2O);
end
% load stats (SNR, FWHM, etc)
load([dataExportPathRef{indexSubj}, 'stats' truncSuffix '.mat'], 'stats');
......@@ -355,12 +449,16 @@ for indexSubj = 1:numberOfSubjects
currentConcentrationRef = concentrationsRef{indexTE,1, 1};
if fitLCModel
% MC_H2O_ConcentrationsRef = concentrations_MC_H2ORef{indexTE,1, 1};
% MC_H2O_ConcentrationsAllSubjRef(indexSubj) = MC_H2O_ConcentrationsRef{2,2};
metaboliteConcentrationsRef = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesLCModel, [], []);
metaboliteConcentrationsRef(end) = metaboliteConcentrationsRef(end) .* 1e-8;
else
metaboliteConcentrationsRef_notScaled = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesProFit, proFitIteration, activeMetabolitesRef);
%scale the concentations similarly as LCModel
metaboliteConcentrationsRef = metaboliteConcentrationsRef_notScaled * 40873 / concentrations_H2ORef * 2;
% MC_H2O_ConcentrationsRef = concentrations_MC_H2ORef * 40873 / concentrations_H2ORef * 2;
% MC_H2O_ConcentrationsAllSubjRef(indexSubj) = MC_H2O_ConcentrationsRef;
end
if ~fitLCModel
currentFQNRef = FQNsRef{indexTE,1, 1};
......@@ -390,17 +488,23 @@ for indexSubj = 1:numberOfSubjects
figure(indexFigure);
hold on
if fitLCModel
% MC_H2O_Concentrations = concentrations_MC_H2O{indexTE,indexAverageDelete, indexCoilDelete};
% MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) = MC_H2O_Concentrations{2,2};
metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesLCModel, [], []);
metaboliteConcentrations(end) = metaboliteConcentrations(end) .* 1e-8;
else
metaboliteConcentrations_notScaled = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesProFit, proFitIteration, activeMetabolites);
%scale the concentations similarly as LCModel
metaboliteConcentrations = metaboliteConcentrations_notScaled * 40873 / concentrations_H2O(runningIndex) * 2;
% MC_H2O_Concentrations = concentrations_MC_H2O(runningIndex) * 40873 / concentrations_H2O(runningIndex) * 2;
% MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) = MC_H2O_Concentrations;
end
if runningIndex <= 3 %scale with coresponding averages
metaboliteConcentrations = metaboliteConcentrations * 3;
MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) = MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) * 3;
else
metaboliteConcentrations = metaboliteConcentrations *3/2;
MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) = MC_H2O_ConcentrationsAllSubj(indexSubj, runningIndex) *3/2;
end
if ~fitLCModel
currentFQN = FQNs{indexTE,indexAverageDelete, indexCoilDelete};
......@@ -448,22 +552,43 @@ for indexSubj = 1:numberOfSubjects
mainMetIterator = iterator(indecesMainMet);
title(SoftwareName)
% if fitLCModel
% title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
% ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
% ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
% num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'...
% ]);
% else
% title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
% ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
% ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
% num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'...
% ', FQN_{Ref}:' num2str(FQN_Ref/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN/max(iterator),3)...
% ', FQN_{Ref}:' num2str(FQN_Ref2/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN2/max(iterator),3)...
% ]);
% end
if fitLCModel
title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'...
]);
else
title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'...
', FQN_{Ref}:' num2str(FQN_Ref/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN/max(iterator),3)...
', FQN_{Ref}:' num2str(FQN_Ref2/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN2/max(iterator),3)...
]);
end
end
end
end
function rpc = calculate_rpc_percent(dataS1, dataS2)
data.set1 = dataS1;
data.set2 = dataS2;
s = size(data.set1);
if ~isequal(s,size(data.set2))
error('data1 and data2 must have the same size');
end
data.set1 = reshape(data.set1, [numel(data.set1),1]);
data.set2 = reshape(data.set2, [numel(data.set2),1]);
data.mask = isfinite(data.set1) & isnumeric(data.set1) & isfinite(data.set2) & isnumeric(data.set2);
data.maskedSet1 = data.set1(data.mask);
data.maskedSet2 = data.set2(data.mask);
data.maskedBaRefData = mean([data.maskedSet1,data.maskedSet2],2);
data.maskedDifferences = (data.maskedSet2-data.maskedSet1) ./ data.maskedBaRefData*100;
differenceSTD = std(data.maskedDifferences, 'omitnan');
rpc = 1.96*differenceSTD;
end
\ No newline at end of file
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment