Plotting error bars on curve (2024)

clear all

%File location of Results to analyze

load ('C:\Users\lexil\Documents\Patient_Study\Healthy Participants\Participant_data/T07b.mat')%young

%If subjects need to be excluded post hoc, you can list the number here

omits = [];

alphaVicon_GC=[];

alpha_IMU=[];

Vicon_GC=[];

IMU_GC=[];

minV = [];

minIMU = [];

maxV = [];

maxIMU = [];

ROMV = [];

ROMRD = [];

RMSE_IMU = [];

p = [];

diff = [];

std_diff = [];

SE = [];

% dat_store{c,1} = gc; %Gait cycle 0-100

% dat_store{c,2} = alphaVicon_GC; %Vicon

% dat_store{c,3} = alpha1D_GC; %RD

% dat_store{c,4} = alpha_MS_GC;

% dat_store{c,5} = ROM;

% dat_store{c,6} = RMSE;

j=0;

sz = size(dat_store);

for i = 1:sz(1)-length(omits) %

if i == omits

else

ROM_all = dat_store{i,4};

ROMIMU_diff(i,:) = mean(ROM_all(:,1:3))-mean(ROM_all(:,4:6));

ROMIMU(i,:) = mean(ROM_all(:,4:6));

ROMVic(i,:) = mean(ROM_all(:,1:3));

alphaVicon_GC = [alphaVicon_GC;(dat_store{i,2})];

%alpha1D_GC = [alpha1D_GC;dat_store{i,3}];

alpha_IMU = [alpha_IMU;dat_store{i,3}];

% Average gc for each individual

Vicon_GC = [Vicon_GC;mean(dat_store{i,2})];

IMU_GC = [IMU_GC;mean(dat_store{i,3})];

RMSE = dat_store{i,5};

%RMSE_1D = [RMSE_1D ;RMSE(2,:)];

RMSE_IMU = [RMSE_IMU ;RMSE(1,:)];

%CORR_1D(i) = corr(reshape(dat_store{i,2},[],1),reshape(dat_store{i,3},[],1));

CORR_IMU(i) = corr(reshape(dat_store{i,2},[],1),reshape(dat_store{i,3},[],1));

end

end

Vicon_GC = Vicon_GC';

IMU_GC = IMU_GC';

diff = Vicon_GC - IMU_GC;

ave_diff = mean(Vicon_GC') - mean(IMU_GC');

std_diff = std(diff,0,2);

SE = std_diff/sqrt(i);

var = ["alphaIMU"];

RMSEtable = [mean(mean(RMSE_IMU,2)) std(mean(RMSE_IMU,2))];

table(var,RMSEtable)

var2 = ["alpha"];

CORR = [mean(CORR_IMU) std(CORR_IMU)];

table(var2,CORR)

GC_all_IMU = reshape(alphaVicon_GC,[],1)-reshape(alpha_IMU,[],1);

GC_all_IMU(abs(GC_all_IMU)<30);

RMSE_IMUa = sqrt(mean(GC_all_IMU.^2))

%ROM

minV = min(ROM_all(:,1));

minIMU = min(ROM_all(:,4));

maxV = max(ROM_all(:,2));

maxIMU = max(ROM_all(:,5));

minstd_V = std(ROM_all(:,1));

minstd_IMU = std(ROM_all(:,4));

maxstd_V = std(ROM_all(:,2));

maxstd_IMU = std(ROM_all(:,5));

ROMstd_V = std(ROM_all(:,3));

ROMstd_IMU = std(ROM_all(:,6));

ROMV = maxV - minV;

ROMRD = maxIMU - minIMU;

gc = dat_store{1,1};

sig = 1.96; %95% Confidence intervals but if sd then replace with 1

figure

patch([gc,flip(gc)],[mean(alphaVicon_GC)-sig*std(alphaVicon_GC) flip(mean(alphaVicon_GC)+sig*std(alphaVicon_GC)) ],[0 0 0],'facealpha',0.2,'edgealpha',0)

hold on

patch([gc,flip(gc)],[mean(alpha_IMU)-sig*std(alpha_IMU) flip(mean(alpha_IMU)+sig*std(alpha_IMU)) ],[1 0 0],'facealpha',0.1,'edgealpha',0)

%patch([gc,flip(gc)],[mean(alpha1D_GC)-sig*std(alpha1D_GC) flip(mean(alpha1D_GC)+sig*std(alpha1D_GC)) ],[0 0 1],'facealpha',0.1,'edgealpha',0)

plot(gc ,mean(alpha_IMU),'r--','LineWidth',1)

%plot(gc ,mean(alpha1D_GC),'b:','LineWidth',1)

plot(gc, mean(alphaVicon_GC),'k-','LineWidth',1)

%xlim([0 TimeEnd])

ylim([-20 80])

grid on

ylabel('Knee flexion angle (\circ)')

xlabel('Gait Cycle (%)')

xticklabels(["0%","20%","40%","60%","80%","100%"])

% g(1) = patch(NaN,NaN,[1 0 0],'facealpha',0.2);

% g(2) = patch(NaN,NaN,[0 0 0],'facealpha',0.2);

% g(3) = patch(NaN,NaN,[0 0 1],'facealpha',0.2);

g(1) = plot(NaN,NaN,'k-','LineWidth',1);

g(2) = plot(NaN,NaN,'r--','LineWidth',1);

rgb = [0 0 0];

FaceAlpha = (0.1);

g(3) = patch([NaN],[NaN],rgb,'EdgeAlpha', 0, 'FaceAlpha',FaceAlpha);

legend(g,'Camera-Marker','IMU','+/- 1.96*SD','location','northwest')

figure()

subplot(3,1,1)

plot(gc ,mean(alpha_IMU),'k--','LineWidth',1)

hold on

plot(gc, mean(alphaVicon_GC),'k-','LineWidth',1)

legend('IMU','Vicon')

xlabel('Gait cycle %')

ylabel('Knee Angle (deg)')

ylim([-20 80]);

grid on

subplot(3,1,2)

plot(gc, (mean(alphaVicon_GC)-mean(alpha_IMU)),'k-','LineWidth',1)

legend('Signed difference')

xlabel('Gait cycle %')

ylabel('Difference (deg)')

ylim([-10 10]);

grid on

subplot(3,1,3)

plot(gc, abs(mean(alphaVicon_GC)-mean(alpha_IMU)),'k-','LineWidth',1)

legend('Absolute difference')

xlabel('Gait cycle %')

ylabel('Difference (deg)')

ylim([-10 10]);

grid on

figure()

x = smooth(resample((gc),101,length(gc),'Dimension',1));

ave_diff_resample = smooth(resample((ave_diff),101,length(ave_diff),'Dimension',1));

SE_resampled = smooth(resample((SE),101,length(SE),'Dimension',1));

set(gcf,'Color','w');

%standard error plot

plot(x, ave_diff_resample,'-','Color',[0 0 0],'linewidth',0.5)

hold on

% Adding error bars

errorbar(x, ave_diff_resample, SE_resampled, 'vertical');

hold on

hline = refline(0, 0);

hline.Color = 'k';

xlim([0 100]);

ylim([-15 15]);

legend('Signed Difference', 'Standard error bar')

xlabel('Gait cycle %')

ylabel('Difference (deg)')

%% Average Bland-Altman plot

% Plot mean difference line

figure()

set(gcf,'Color','w');

ave_vic = mean(Vicon_GC');

ave_imu = mean(IMU_GC');

mean_data = (ave_vic + ave_imu)/2;

mean_data = resample((mean_data),101,length(mean_data),'Dimension',1);

std_diff = resample((std_diff),101,length(std_diff),'Dimension',1);

first_index = 1;

last_index = numel(mean_data);

scatter(mean_data(first_index), ave_diff_resample(first_index),80, '*');

hold on

scatter(mean_data(last_index), ave_diff_resample(last_index),80, '^');

hold on

scatter(mean_data(2:100), ave_diff_resample(2:100), 'filled','MarkerFaceColor', 'black');

hold on;

%+- 2SD

plot([0 120], [mean(ave_diff_resample)'+ 2*mean(std_diff), mean(ave_diff_resample)' + 2*mean(std_diff)], 'r--');

plot([0 120], [mean(ave_diff_resample)' - 2*mean(std_diff), mean(ave_diff_resample)' - 2*mean(std_diff)], 'r--');

%mean diff line

plot([0 120], [mean(ave_diff_resample), mean(ave_diff_resample)], 'k-', 'LineWidth', 1);

%calculate error bars

%error bars are the standard error of the sample

error_bar_y = errorbar(mean_data,ave_diff_resample,SE_resampled, 'vertical', 'LineStyle', 'none');

set(error_bar_y, 'color', 'k', 'LineWidth', 0.5)

legend({'HS', 'HS', 'Average gait cycle Data', '+ 2 SD', '- 2 SD', 'Mean Difference', 'SE'});

hold on

xlim([0 120])

ylim([-10 10])

xlabel('Mean($\o - \bar{\o}$) (deg)', 'Interpreter','Latex')

ylabel('Difference (deg)')

title('Cycling - Younger Adults')