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

Changed a lot of things in plot_all_test_data_profit

parent e401d52a
No related branches found
No related tags found
No related merge requests found
function plot_all_test_data_profit() function plot_all_test_data_profit()
fitLCModel = true; %false fits profit
reconOption = 5; reconOption = 5;
truncOn = false; truncOn = false;
...@@ -13,7 +12,7 @@ truncOn = false; ...@@ -13,7 +12,7 @@ truncOn = false;
offsetScatter = 0; offsetScatter = 0;
markerScatter = 'd'; markerScatter = 'd';
markerColor = [0 1 0]; markerColor = [0 1 0];
% %
reconOption = 2; reconOption = 2;
truncOn = true; truncOn = true;
offsetScatter = 0.15; offsetScatter = 0.15;
...@@ -21,16 +20,6 @@ markerScatter = 'x'; ...@@ -21,16 +20,6 @@ markerScatter = 'x';
markerColor = [1 0 0]; markerColor = [1 0 0];
if fitLCModel
SoftwareName = 'LCModel';
indexFigSoftware = 0;
else
SoftwareName = 'ProFit';
indexFigSoftware = 50;
proFitIteration =4;
end
plotInSameFigureAllSubjects = true;
exportFolders = {'0 Water References',... exportFolders = {'0 Water References',...
'1 Normal recon', '2 Removed averages ','3 Removed coils',... '1 Normal recon', '2 Removed averages ','3 Removed coils',...
'4 Removed averages and coils', '5 Water Supp Off', '6 Frequency Alignment Off', ... '4 Removed averages and coils', '5 Water Supp Off', '6 Frequency Alignment Off', ...
...@@ -38,8 +27,8 @@ exportFolders = {'0 Water References',... ...@@ -38,8 +27,8 @@ exportFolders = {'0 Water References',...
pathNameExport = 'ProFit test data'; pathNameExport = 'ProFit test data';
pathBaseExportFiles = pathToDataFolder(pathNameExport, exportFolders); pathBaseExportFiles = pathToDataFolder(pathNameExport, exportFolders);
pathBaseExportFiles = pathBaseExportFiles(1:end-1); % 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 trunc\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 sinebell\']; % pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 sinebell\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 matched\']; % pathBaseExportFiles = [pathBaseExportFiles ' - Weighted Cost Fun R3 0.25scale m_15 matched\'];
% pathBaseExportFiles = [pathBaseExportFiles ' - R1 0.25scale m_15 trunc\']; % pathBaseExportFiles = [pathBaseExportFiles ' - R1 0.25scale m_15 trunc\'];
...@@ -110,6 +99,194 @@ for indexMet = 1:numberOfMet ...@@ -110,6 +99,194 @@ for indexMet = 1:numberOfMet
end end
indecesMainMet = logical(indecesMainMet); indecesMainMet = logical(indecesMainMet);
fitLCModel = true; %false fits profit
if fitLCModel
SoftwareName = 'LCModel';
else
SoftwareName = 'ProFit';
end
[metConcAllSubjLCModel, metConcAllSubjRefLCModel, metConcAllSubjNormLCModel] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
offsetScatter, markerScatter, markerColor);
fitLCModel = false; %false fits profit
if fitLCModel
SoftwareName = 'LCModel';
else
SoftwareName = 'ProFit';
end
[metConcAllSubjProFit, metConcAllSubjRefProFit, metConcAllSubjNormProFit] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
offsetScatter, markerScatter, markerColor);
metaboliteRef = 14;%'tCho+' 3;%'tCr(CH3)' 12;%'NAA(CH3)'
metConcAllSubjProFitScaled = metConcAllSubjProFit ./ metConcAllSubjProFit(:,:,metaboliteRef);
metConcAllSubjRefProFitScaled = metConcAllSubjRefProFit ./ metConcAllSubjRefProFit(:,metaboliteRef);
metConcAllSubjLCModelScaled = metConcAllSubjLCModel ./ metConcAllSubjLCModel(:,:,metaboliteRef);
metConcAllSubjRefLCModelScaled = metConcAllSubjRefLCModel ./ metConcAllSubjRefLCModel(:,metaboliteRef);
metaboliteRef = 0;
metConcAllSubjProFitScaled = metConcAllSubjProFit;
metConcAllSubjRefProFitScaled = metConcAllSubjRefProFit;
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'};
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_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
limits = 'auto'; % how to set the axes limits
if 1 % colors for the data sets may be set as:
colors = 'br'; % character codes
else
colors = [0 0 1;... % or RGB triplets
1 0 0];
end
numOfSimulations = length(namingSuffixesAve)*length(namingSuffixesCoil);
rpc_table = cell(numberOfMet+3,3);
rpc_table{1,2} = softwareNames{1};
rpc_table{1,3} = softwareNames{2};
rpc_table(2:end-2,1) = metabolitesToDisplay;
rpc_table{end,1} = 'Mean';
for indexMet = 1: numberOfMet
if indexMet == metaboliteRef
continue;
end
%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
%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
%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,...
% '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,...
% '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,...
% '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;
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;
xlswrite([pathBaseExportFiles, 'InVivoResults_RPC.xlsx'], rpc_table)
meanMetDeviation = mean(abs(metaboliteConcentrationsAllSubjNorm), 'omitnan');
stdMetDeviation = std(abs(metaboliteConcentrationsAllSubjNorm), 'omitnan');
meanAllMetDeviation = mean(abs(metaboliteConcentrationsAllSubjNorm(:)), 'omitnan');
stdAllMetDeviation = std(abs(metaboliteConcentrationsAllSubjNorm(:)), 'omitnan');
mainMetaboliteConcentrationsAllSubj = metaboliteConcentrationsAllSubjNorm(:, indecesMainMet);
meanMainMetDeviation = mean(abs(mainMetaboliteConcentrationsAllSubj(:)), 'omitnan');
stdMainMetDeviation = std(abs(mainMetaboliteConcentrationsAllSubj(:)), 'omitnan');
metaboliteTable = cell(numberOfMet+3,2);
metaboliteTable(2:numberOfMet+1,1) = metabolitesLabels;
metaboliteTable{numberOfMet+2,1} = 'Mean';
metaboliteTable{numberOfMet+3,1} = 'Mean Main Metabolites';
plusMinusSign = char(177);
metaboliteTable{1,2} = SoftwareName;
for indexMet = 1:numberOfMet
metaboliteTable{indexMet+1,2} = [num2str(meanMetDeviation(indexMet),'%.1f'), plusMinusSign, num2str(stdMetDeviation(indexMet),'%.1f')];
end
metaboliteTable{numberOfMet+2,2} = [num2str(meanAllMetDeviation,'%.1f'), plusMinusSign, num2str(stdAllMetDeviation,'%.1f')];
metaboliteTable{numberOfMet+3,2} = [num2str(meanMainMetDeviation,'%.1f'), plusMinusSign, num2str(stdMainMetDeviation,'%.1f')];
xlswrite([pathBaseExportFiles, 'InVivoResults' SoftwareName, '.xlsx'], metaboliteTable)
%
% mean(cell2mat(NAA_SNR_all(:,1)))
% std(cell2mat(NAA_SNR_all(:,1)))
% mean(cell2mat(NAA_SNR_all(:,2)))
% mean(cell2mat(NAA_SNR_all(:,3)))
% mean(cell2mat(NAA_SNR_all(:,4)))
% mean(cell2mat(NAA_SNR_all(:,5)))
% std(cell2mat(NAA_SNR_all(:,5)))
end
function metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesNames, proFitIteration, activeMetabolites)
metaboliteConcentrations = zeros(size(metabolitesNames));
if fitLCModel
activeMetabolites = currentConcentration(1,:);
currentConcentration = currentConcentration(2,:);
else
activeMetabolites = activeMetabolites{proFitIteration};
currentConcentration = currentConcentration{proFitIteration};
end
for indexMet = 1 : length(metabolitesNames)
indexMetFit = find(strcmp(activeMetabolites,metabolitesNames(indexMet)));
if fitLCModel
metaboliteConcentrations(indexMet) = currentConcentration{indexMetFit};
else
metaboliteConcentrations(indexMet) = currentConcentration(indexMetFit);
end
end
% metaboliteConcentrations(1) = metaboliteConcentrations(1) / 1e6; % adjust the MMB to be within range
end
function [metaboliteConcentrationsAllSubj, metaboliteConcentrationsAllSubjRef, metaboliteConcentrationsAllSubjNorm] = ...
getAllMetaboliteConc(dataExportPathBase, dataExportPathBaseRef, truncSuffix, subjects, TEs, ...
namingSuffixesAve, namingSuffixesCoil, fitLCModel, reconOption, ...
numberOfMet, metabolitesLCModel, metabolitesProFit, metabolitesLabels, indecesMainMet, ...
offsetScatter, markerScatter, markerColor)
% settings
plotInSameFigureAllSubjects = true;
scalingFactorColor = 1;
yLimValue = 50;
NAA_SNR_global_ref = 700;% 650; for truncOn true
if fitLCModel
SoftwareName = 'LCModel';
indexFigSoftware = 0;
else
SoftwareName = 'ProFit';
indexFigSoftware = 50;
proFitIteration =4;
end
% prepare arrays etc
numberOfSubjects = length(subjects); numberOfSubjects = length(subjects);
dataExportPathRef = cell(1,numberOfSubjects); dataExportPathRef = cell(1,numberOfSubjects);
dataExportPath = cell(1,numberOfSubjects); dataExportPath = cell(1,numberOfSubjects);
...@@ -117,18 +294,12 @@ for indexSubj = 1:numberOfSubjects ...@@ -117,18 +294,12 @@ for indexSubj = 1:numberOfSubjects
dataExportPathRef{indexSubj} = [dataExportPathBaseRef subjects{indexSubj} '\\']; dataExportPathRef{indexSubj} = [dataExportPathBaseRef subjects{indexSubj} '\\'];
dataExportPath{indexSubj} = [dataExportPathBase subjects{indexSubj} '\\']; dataExportPath{indexSubj} = [dataExportPathBase subjects{indexSubj} '\\'];
end end
numberOfTEs = length(TEs);
figId = figure(); figId = figure();
colormapVect = colormap(); colormapVect = colormap();
close(figId); close(figId);
colormapVect = colormapVect(end:-1:1, :); colormapVect = colormapVect(end:-1:1, :);
numberOfTEs = length(TEs);
numberOfColors = size(colormapVect,1); numberOfColors = size(colormapVect,1);
scalingFactorColor = 1;
yLimValue = 50;
NAA_SNR_global_ref = 700;% 650; for truncOn true
%% do the actual plotting
NAA_SNR_all = cell(numberOfSubjects, numberOfTEs); NAA_SNR_all = cell(numberOfSubjects, numberOfTEs);
meanPercentualChange = 0; meanPercentualChange = 0;
iterator = 0; iterator = 0;
...@@ -136,7 +307,10 @@ nansIterator = 0; ...@@ -136,7 +307,10 @@ nansIterator = 0;
meanFQN = 0; meanFQN = 0;
meanFQN2 = 0; meanFQN2 = 0;
numOfSimulations = length(namingSuffixesAve)*length(namingSuffixesCoil); numOfSimulations = length(namingSuffixesAve)*length(namingSuffixesCoil);
metaboliteConcentrationsAllSubj = zeros(numberOfSubjects * numOfSimulations, numberOfMet); metaboliteConcentrationsAllSubjNorm = zeros(numberOfSubjects * numOfSimulations, numberOfMet);
metaboliteConcentrationsAllSubj = zeros(numberOfSubjects, numOfSimulations, numberOfMet);
metaboliteConcentrationsAllSubjRef = zeros(numberOfSubjects, numberOfMet);
%% do the actual plotting
for indexSubj = 1:numberOfSubjects for indexSubj = 1:numberOfSubjects
if fitLCModel if fitLCModel
load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations'); load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations');
...@@ -144,11 +318,15 @@ for indexSubj = 1:numberOfSubjects ...@@ -144,11 +318,15 @@ for indexSubj = 1:numberOfSubjects
load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations'); load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '.mat'], 'concentrations');
else else
load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2'); load([dataExportPathRef{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2');
concentrationsRef = concentrations; load([dataExportPathRef{indexSubj}, 'concentrations', truncSuffix '_H2O_profit.mat'], 'concentrations_H2O');
concentrationsRef = concentrations;
concentrations_H2ORef = cell2mat(concentrations_H2O);
activeMetabolitesRef = activeMetabolites; activeMetabolitesRef = activeMetabolites;
FQNsRef = FQNs; FQNsRef = FQNs;
FQNsRef2 = FQNs2; FQNsRef2 = FQNs2;
load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2'); load([dataExportPath{indexSubj}, 'concentrations' truncSuffix '_profit.mat'], 'concentrations', 'activeMetabolites', 'FQNs', 'FQNs2');
load([dataExportPath{indexSubj}, 'concentrations', truncSuffix '_H2O_profit.mat'], 'concentrations_H2O');
concentrations_H2O = cell2mat(concentrations_H2O);
end end
% load stats (SNR, FWHM, etc) % load stats (SNR, FWHM, etc)
load([dataExportPathRef{indexSubj}, 'stats' truncSuffix '.mat'], 'stats'); load([dataExportPathRef{indexSubj}, 'stats' truncSuffix '.mat'], 'stats');
...@@ -178,8 +356,11 @@ for indexSubj = 1:numberOfSubjects ...@@ -178,8 +356,11 @@ for indexSubj = 1:numberOfSubjects
if fitLCModel if fitLCModel
metaboliteConcentrationsRef = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesLCModel, [], []); metaboliteConcentrationsRef = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesLCModel, [], []);
metaboliteConcentrationsRef(end) = metaboliteConcentrationsRef(end) .* 1e-8;
else else
metaboliteConcentrationsRef = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesProFit, proFitIteration, activeMetabolitesRef); metaboliteConcentrationsRef_notScaled = matchFittedMetabolites(fitLCModel, currentConcentrationRef, metabolitesProFit, proFitIteration, activeMetabolitesRef);
%scale the concentations similarly as LCModel
metaboliteConcentrationsRef = metaboliteConcentrationsRef_notScaled * 40873 / concentrations_H2ORef * 2;
end end
if ~fitLCModel if ~fitLCModel
currentFQNRef = FQNsRef{indexTE,1, 1}; currentFQNRef = FQNsRef{indexTE,1, 1};
...@@ -188,14 +369,14 @@ for indexSubj = 1:numberOfSubjects ...@@ -188,14 +369,14 @@ for indexSubj = 1:numberOfSubjects
FQN_Ref2 = currentFQNRef2{proFitIteration}; FQN_Ref2 = currentFQNRef2{proFitIteration};
end end
%scale with the MC water reference amplitude %scale with the MC water reference amplitude
metaboliteConcentrationsRef = metaboliteConcentrationsRef ./ water_ref_amplitude;
metaboliteConcentrationsRefNormalized = (metaboliteConcentrationsRef-metaboliteConcentrationsRef) ./ metaboliteConcentrationsRef *100; metaboliteConcentrationsRefNormalized = (metaboliteConcentrationsRef-metaboliteConcentrationsRef) ./ metaboliteConcentrationsRef *100;
scatter([1:numberOfMet]-offsetScatter,metaboliteConcentrationsRefNormalized,[],'d', 'MarkerEdgeColor', currentColorRef); scatter([1:numberOfMet]-offsetScatter,metaboliteConcentrationsRefNormalized,[],'d', 'MarkerEdgeColor', currentColorRef);
% cb = colorbar; metaboliteConcentrationsAllSubjRef(indexSubj, :) = metaboliteConcentrationsRef;
% cb.Limits = [0 NAA_ref_SNR*scalingFactorColor]; % cb = colorbar;
% cb.Color = colormapVect % cb.Limits = [0 NAA_ref_SNR*scalingFactorColor];
% cb.Color = colormapVect
metaboliteConcentrationsAll = zeros(numOfSimulations, numberOfMet); metaboliteConcentrationsAll = zeros(numOfSimulations, numberOfMet);
% do plotting for each individual case % do plotting for each individual case
for indexAverageDelete = 1: length(namingSuffixesAve) for indexAverageDelete = 1: length(namingSuffixesAve)
...@@ -210,8 +391,16 @@ for indexSubj = 1:numberOfSubjects ...@@ -210,8 +391,16 @@ for indexSubj = 1:numberOfSubjects
hold on hold on
if fitLCModel if fitLCModel
metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesLCModel, [], []); metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesLCModel, [], []);
else metaboliteConcentrations(end) = metaboliteConcentrations(end) .* 1e-8;
metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesProFit, proFitIteration, activeMetabolites); else
metaboliteConcentrations_notScaled = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesProFit, proFitIteration, activeMetabolites);
%scale the concentations similarly as LCModel
metaboliteConcentrations = metaboliteConcentrations_notScaled * 40873 / concentrations_H2O(runningIndex) * 2;
end
if runningIndex <= 3 %scale with coresponding averages
metaboliteConcentrations = metaboliteConcentrations * 3;
else
metaboliteConcentrations = metaboliteConcentrations *3/2;
end end
if ~fitLCModel if ~fitLCModel
currentFQN = FQNs{indexTE,indexAverageDelete, indexCoilDelete}; currentFQN = FQNs{indexTE,indexAverageDelete, indexCoilDelete};
...@@ -220,26 +409,24 @@ for indexSubj = 1:numberOfSubjects ...@@ -220,26 +409,24 @@ for indexSubj = 1:numberOfSubjects
currentFQN2 = FQNs2{indexTE,indexAverageDelete, indexCoilDelete}; currentFQN2 = FQNs2{indexTE,indexAverageDelete, indexCoilDelete};
FQN2 = currentFQN2{proFitIteration}; FQN2 = currentFQN2{proFitIteration};
meanFQN2 = meanFQN2 + FQN2; meanFQN2 = meanFQN2 + FQN2;
end end
%scale with the MC water reference amplitude
metaboliteConcentrations = metaboliteConcentrations ./ water_amplitude;
metaboliteConcentrationsNormalized = (metaboliteConcentrations-metaboliteConcentrationsRef) ... metaboliteConcentrationsNormalized = (metaboliteConcentrations-metaboliteConcentrationsRef) ...
./ metaboliteConcentrationsRef *100; ./ metaboliteConcentrationsRef *100;
metaboliteConcentrationsAll(runningIndex, :) = metaboliteConcentrationsNormalized; metaboliteConcentrationsAll(runningIndex, :) = metaboliteConcentrationsNormalized;
metaboliteConcentrationsNormalized(abs(metaboliteConcentrationsNormalized) > yLimValue*2) = NaN; %TODO change! This is actually stupid metaboliteConcentrationsNormalized(abs(metaboliteConcentrationsNormalized) > yLimValue*2) = NaN; %TODO change! This is actually stupid
metaboliteConcentrationsAllSubj(runningIndex + (indexSubj-1)*numOfSimulations, :) = metaboliteConcentrationsNormalized; metaboliteConcentrationsAllSubjNorm(runningIndex + (indexSubj-1)*numOfSimulations, :) = metaboliteConcentrationsNormalized;
metaboliteConcentrationsAllSubj(indexSubj, runningIndex, :) = metaboliteConcentrations;
if reconOption>4 if reconOption>4
%use this scatter if you want to check if there are structural offsets %use this scatter if you want to check if there are structural offsets
scatter([1:numberOfMet]+offsetScatter,metaboliteConcentrationsNormalized, markerScatter, 'MarkerEdgeColor', markerColor, 'LineWidth', 1); % scatter([1:numberOfMet]+offsetScatter,metaboliteConcentrationsNormalized, markerScatter, 'MarkerEdgeColor', markerColor, 'LineWidth', 1);
else else
% SNR matched % SNR matched
currentColor = colormapVect(round(NAA_SNR/NAA_SNR_global_ref/scalingFactorColor * numberOfColors),:); currentColor = colormapVect(round(NAA_SNR/NAA_SNR_global_ref/scalingFactorColor * numberOfColors),:);
scatter([1:numberOfMet]+offsetScatter,metaboliteConcentrationsNormalized, markerScatter, 'MarkerEdgeColor',currentColor, 'LineWidth', 2); % scatter([1:numberOfMet]+offsetScatter,metaboliteConcentrationsNormalized, markerScatter, 'MarkerEdgeColor',currentColor, 'LineWidth', 2);
end end
currentNans = isnan(metaboliteConcentrationsNormalized); currentNans = isnan(metaboliteConcentrationsNormalized);
nansIterator = nansIterator + currentNans; nansIterator = nansIterator + currentNans;
meanPercentualChange = meanPercentualChange + abs(metaboliteConcentrationsNormalized); meanPercentualChange = meanPercentualChange + abs(metaboliteConcentrationsNormalized);
iterator = iterator + ~(currentNans); iterator = iterator + ~(currentNans);
end end
end end
...@@ -248,92 +435,35 @@ for indexSubj = 1:numberOfSubjects ...@@ -248,92 +435,35 @@ for indexSubj = 1:numberOfSubjects
cb = colorbar('XTickLabel',num2str(round(xticksSNR)'), 'Direction','reverse'); cb = colorbar('XTickLabel',num2str(round(xticksSNR)'), 'Direction','reverse');
ylabel(cb,'SNR_{NAA(CH3)}'); ylabel(cb,'SNR_{NAA(CH3)}');
% boxplot(metaboliteConcentrationsAll); % boxplot(metaboliteConcentrationsAll);
xticks(1:numberOfMet); xticks(1:numberOfMet);
xticklabels(metabolitesLabels); xticklabels(metabolitesLabels);
xtickangle(45) xtickangle(45)
ylim([-yLimValue,yLimValue]); ylim([-yLimValue,yLimValue]);
xlim([0.5,17.5]) xlim([0.5,17.5])
ylabel('Concentration change c_{k,%} (%)') ylabel('Concentration change c_{k,%} (%)')
mainMetMeanChange = meanPercentualChange(indecesMainMet); mainMetMeanChange = meanPercentualChange(indecesMainMet);
mainMetNaNs = nansIterator(indecesMainMet); mainMetNaNs = nansIterator(indecesMainMet);
mainMetIterator = iterator(indecesMainMet); mainMetIterator = iterator(indecesMainMet);
title(SoftwareName) title(SoftwareName)
% if fitLCModel % if fitLCModel
% title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'... % title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
% ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '... % ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
% ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '... % ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
% num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'... % num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'...
% ]); % ]);
% else % else
% title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'... % title([SoftwareName ... ' Subject: ', subjects{indexSubj}, ' TE ', num2str(TEs(indexTE)) ' ms'...
% ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '... % ' - mean change (main mets): ', num2str(mean(meanPercentualChange./iterator,'omitnan'),3) '%; '...
% ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '... % ' (', num2str(mean(mainMetMeanChange./mainMetIterator,'omitnan'),3) '%); '...
% num2str(sum(nansIterator)), ' (', num2str(sum(mainMetNaNs)), ') NaNs'... % 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_Ref/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN/max(iterator),3)...
% ', FQN_{Ref}:' num2str(FQN_Ref2/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN2/max(iterator),3)... % ', FQN_{Ref}:' num2str(FQN_Ref2/numberOfSubjects,3), ', FQN_{mean}:' num2str(meanFQN2/max(iterator),3)...
% ]); % ]);
% end % end
end
end
meanMetDeviation = mean(abs(metaboliteConcentrationsAllSubj), 'omitnan');
stdMetDeviation = std(abs(metaboliteConcentrationsAllSubj), 'omitnan');
meanAllMetDeviation = mean(abs(metaboliteConcentrationsAllSubj(:)), 'omitnan');
stdAllMetDeviation = std(abs(metaboliteConcentrationsAllSubj(:)), 'omitnan');
mainMetaboliteConcentrationsAllSubj = metaboliteConcentrationsAllSubj(:, indecesMainMet);
meanMainMetDeviation = mean(abs(mainMetaboliteConcentrationsAllSubj(:)), 'omitnan');
stdMainMetDeviation = std(abs(mainMetaboliteConcentrationsAllSubj(:)), 'omitnan');
metaboliteTable = cell(numberOfMet+3,2);
metaboliteTable(2:numberOfMet+1,1) = metabolitesLabels;
metaboliteTable{numberOfMet+2,1} = 'Mean';
metaboliteTable{numberOfMet+3,1} = 'Mean Main Metabolites';
plusMinusSign = char(177);
metaboliteTable{1,2} = SoftwareName;
for indexMet = 1:numberOfMet
metaboliteTable{indexMet+1,2} = [num2str(meanMetDeviation(indexMet),'%.1f'), plusMinusSign, num2str(stdMetDeviation(indexMet),'%.1f')];
end
metaboliteTable{numberOfMet+2,2} = [num2str(meanAllMetDeviation,'%.1f'), plusMinusSign, num2str(stdAllMetDeviation,'%.1f')];
metaboliteTable{numberOfMet+3,2} = [num2str(meanMainMetDeviation,'%.1f'), plusMinusSign, num2str(stdMainMetDeviation,'%.1f')];
xlswrite([pathBaseExportFiles, 'InVivoResults' SoftwareName, '.xlsx'], metaboliteTable)
%
% mean(cell2mat(NAA_SNR_all(:,1)))
% std(cell2mat(NAA_SNR_all(:,1)))
% mean(cell2mat(NAA_SNR_all(:,2)))
% mean(cell2mat(NAA_SNR_all(:,3)))
% mean(cell2mat(NAA_SNR_all(:,4)))
% mean(cell2mat(NAA_SNR_all(:,5)))
% std(cell2mat(NAA_SNR_all(:,5)))
end
function metaboliteConcentrations = matchFittedMetabolites(fitLCModel, currentConcentration, metabolitesNames, proFitIteration, activeMetabolites)
metaboliteConcentrations = zeros(size(metabolitesNames));
if fitLCModel
activeMetabolites = currentConcentration(1,:);
currentConcentration = currentConcentration(2,:);
else
activeMetabolites = activeMetabolites{proFitIteration};
currentConcentration = currentConcentration{proFitIteration};
end
for indexMet = 1 : length(metabolitesNames)
indexMetFit = find(strcmp(activeMetabolites,metabolitesNames(indexMet)));
if fitLCModel
metaboliteConcentrations(indexMet) = currentConcentration{indexMetFit};
else
metaboliteConcentrations(indexMet) = currentConcentration(indexMetFit);
end end
end end
% metaboliteConcentrations(1) = metaboliteConcentrations(1) / 1e6; % adjust the MMB to be within range
end end
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment