Calculations/Data-Analyzer/+Analyzer/conductPCA.m

112 lines
4.9 KiB
Matlab

function conductPCA(od_imgs, scan_reference_values, scan_parameter_values, doPlot, doSave, saveDir)
%% Performs PCA on optical density images, visualizes and optionally saves results.
%
% Inputs:
% od_imgs - cell array of OD images
% scan_reference_values - array of unique control parameter values
% scan_parameter_values - array mapping each image to a control parameter
% doPlot - logical, true to plot figures
% doSave - logical, true to save figures
% saveDir - directory to save figures if doSave is true
%
% Requires:
% +Calculator/computeCumulants.m
if nargin < 4, doPlot = true; end
if nargin < 5, doSave = false; end
if nargin < 6, saveDir = pwd; end
%% PCA computation
allImgs3D = cat(3, od_imgs{:});
[Nx, Ny] = size(allImgs3D(:,:,1));
Xall = reshape(allImgs3D, [], numel(od_imgs))';
[coeff, score, ~, ~, explained] = pca(Xall);
figCount = 1;
%% --- Figure 1: PC1 Image ---
if doPlot
pc1_vector = coeff(:,1);
pc1_image = reshape(pc1_vector, Nx, Ny);
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]);
imagesc(pc1_image); axis image off; colormap(Colormaps.coolwarm()); colorbar;
title(sprintf('First Principal Component (PC1) Image - Explains %.2f%% Variance', explained(1)));
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_Image.png')); end
figCount = figCount + 1;
end
%% --- Figure 2: PC1 Scores Scatter ---
if doPlot
numGroups = numel(scan_reference_values);
colors = lines(numGroups);
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]); hold on;
for g = 1:numGroups
idx = scan_parameter_values == scan_reference_values(g);
scatter(repmat(scan_reference_values(g), sum(idx),1), score(idx,1), 36, colors(g,:), 'filled');
end
xlabel('Control Parameter'); ylabel('PC1 Score');
title('Evolution of PC1 Scores'); grid on;
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_Scatter.png')); end
figCount = figCount + 1;
end
%% --- Figure 3: PC1 Distributions ---
if doPlot
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]);
hold on;
for g = 1:numGroups
idx = scan_parameter_values == scan_reference_values(g);
data = score(idx,1);
histogram(data, 'Normalization', 'pdf', 'FaceColor', colors(g,:), 'FaceAlpha', 0.3);
[f, xi] = ksdensity(data);
plot(xi, f, 'Color', colors(g,:), 'LineWidth', 2);
end
xlabel('PC1 Score'); ylabel('Probability Density');
title('PC1 Score Distributions by Group');
legend(arrayfun(@num2str, scan_reference_values, 'UniformOutput', false), 'Location', 'Best');
grid on;
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_Distributions.png')); end
figCount = figCount + 1;
end
%% --- Figure 4: Boxplot of PC1 Scores ---
if doPlot
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]);
boxplot(score(:,1), scan_parameter_values);
xlabel('Control Parameter'); ylabel('PC1 Score');
title('PC1 Score Boxplots by Group'); grid on;
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_Boxplot.png')); end
figCount = figCount + 1;
end
%% --- Figure 5: Mean ± SEM of PC1 Scores ---
if doPlot
meanScores = arrayfun(@(g) mean(score(scan_parameter_values == g,1)), scan_reference_values);
semScores = arrayfun(@(g) std(score(scan_parameter_values == g,1))/sqrt(sum(scan_parameter_values == g)), scan_reference_values);
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]);
errorbar(scan_reference_values, meanScores, semScores, '-o', 'LineWidth', 2);
xlabel('Control Parameter'); ylabel('Mean PC1 Score ± SEM');
title('Mean ± SEM of PC1 Scores'); grid on;
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_MeanSEM.png')); end
figCount = figCount + 1;
end
%% --- Figure 6: Binder Cumulant ---
if doPlot
binderVals = arrayfun(@(g) ...
Calculator.computeCumulants(score(scan_parameter_values == g,1)), ...
scan_reference_values);
figure(figCount); clf; set(gcf, 'Color', 'w', 'Position', [100 100 950 750]);
plot(scan_reference_values, binderVals, '-o', 'LineWidth', 2);
xlabel('Control Parameter'); ylabel('Binder Cumulant (PC1)');
title('Binder Cumulant of PC1 Scores'); grid on;
if doSave, saveas(gcf, fullfile(saveDir, 'PC1_BinderCumulant.png')); end
end
%% --- ANOVA Test ---
p = anova1(score(:,1), arrayfun(@num2str, scan_parameter_values, 'UniformOutput', false), 'off');
fprintf('[INFO] ANOVA p-value for PC1 score differences between groups: %.4e\n', p);
end