Skip to content

Commit

Permalink
Individual data points for Figure 3,4 and 5
Browse files Browse the repository at this point in the history
  • Loading branch information
chrsi30 committed Sep 6, 2021
1 parent b2d9179 commit 69f73fb
Show file tree
Hide file tree
Showing 17 changed files with 229 additions and 108 deletions.
Binary file removed Bergqvist et al (2017)/GLUTPlosOne.mexw64
Binary file not shown.
Binary file removed Bergqvist et al (2017)/GLUTPlosOnekr.mexw64
Binary file not shown.
Binary file removed Bergqvist et al (2017)/GLUTcomplete.mexw64
Binary file not shown.
Binary file removed Bergqvist et al (2017)/GLUTcompletekr.mexw64
Binary file not shown.
Binary file modified MLMM/Data/HFD_BW_8W.mat
Binary file not shown.
Binary file added MLMM/Data/RAW_BW_DATA.mat
Binary file not shown.
Binary file added MLMM/Data/RAW_IC_DATA.mat
Binary file not shown.
Binary file added MLMM/Data/RAW_IPGTT_DATA.mat
Binary file not shown.
7 changes: 7 additions & 0 deletions MLMM/Data/Setup_data.m
Original file line number Diff line number Diff line change
Expand Up @@ -39,4 +39,11 @@
D.Phos.SEMPKB= MeanSEM_SR2019_PKB_IRS1_AS160.SEM_pPKB;
D.Phos.MeanAS160= MeanSEM_SR2019_PKB_IRS1_AS160.Mean_pAS160;
D.Phos.SEMAS160= MeanSEM_SR2019_PKB_IRS1_AS160.SEM_pAS160;

load('RAW_BW_DATA.mat')
D.RAW_BW_DATA = RAW_BW_DATA;
load('RAW_IPGTT_DATA.mat')
D.RAW_IPGTT_DATA = RAW_IPGTT_DATA;
load('RAW_IC_DATA.mat')
D.RAW_IC_DATA = RAW_IC_DATA;
clear("Stenkula*")
208 changes: 153 additions & 55 deletions MLMM/Figure_3_and_4.m
Original file line number Diff line number Diff line change
Expand Up @@ -3,25 +3,24 @@
id_FFM = ismember(SBstates(objModel),'FFM') ; id_FM = ismember(SBstates(objModel),'FM');
id_FFMname = ismember(pNames,'FFMinit') ; id_FMname = ismember(pNames,'FMinit');
t_diet_up=0:1:56;
t_fasting_up=[ 42 42+0.5833 ] ; %t_diet:0.0001:t_diet+0.5833;
t_fasting_up=[ 42 42+0.5833 ] ;
t_IPGTT_up=(0:1:120)/(60*24)+t_fasting_up(end);

t_diet=14;
t_fasting=[ t_diet t_diet+0.5833 ] ; %t_diet:0.0001:t_diet+0.5833;
t_fasting=[ t_diet t_diet+0.5833 ] ;
t_IPGTT=(0:1:120)/(60*24)+t_fasting(end);

id_BW=ismember(SBvariables(objModel),'BW');
id_Glu = ismember(SBvariables(objModel),'Gconc');
id_Ins = ismember(SBvariables(objModel),'Iconc');

%% UNC
%% Model uncertainty
for i = 1:nfiles
load( strcat( filedir ,'/',files(i).name ) )
disp(strcat('Processing - ',files(i).name, int2str(i),'/',int2str(nfiles)))
theta=Results.xbest;
theta=Results.xbest;


[ AllSimulations ] = Function_AllSimulations( theta ,func_mex_model, D, utility);
[ AllSimulations ] = Function_AllSimulations( theta ,func_mex_model, D, utility); % Model simulation

chow8w = AllSimulations.chow8w;
hfd8w = AllSimulations.hfd8w;
Expand Down Expand Up @@ -154,7 +153,6 @@
end

%% Plotting

predColor = [0.9,.9,0.9];
TrainColor = [ 0.5,0.5,0.5];
label_size = 25;
Expand Down Expand Up @@ -182,17 +180,16 @@
best_cost_pred_data = x_pred ;
disp(strcat('Best cost for validation data = ', num2str(best_cost_pred_data)));



%% Figure 1 - Body weight - TRAINING AND PREDICTIONS
%% Figure 3 - Body weight - TRAINING AND PREDICTIONS
t_diet=0:1:14;
fill_time_2w = [ t_diet fliplr(t_diet) ];

fig=fig+1;
figure (fig)
set(figure(fig), 'outerposition',[0 0 2560 1440], 'PaperType','a4')

%%% Figure 1a BW - Hansson 2018 et al
%Figure 3a-d

Experiments_to_plot = {'chow2w','chow10dHFD4d','chow8dHFD6d','hfd2w'};

for kk = 1:length(Experiments_to_plot)
Expand All @@ -214,25 +211,49 @@
time_plot = t_diet;
time_fill = fill_time_2w;
end


f=fill(time_fill, [ values(:,2)' fliplr( values(:,1)' ) ], TrainColor, 'EdgeColor', 'none' );
set(f,'facealpha',.5)
hold on

if kk == 1
Raw_data1 = D.RAW_BW_DATA.chow_2w_basal ; Raw_data2 = D.RAW_BW_DATA.chow_2w_end;
elseif kk == 2
Raw_data1 = D.RAW_BW_DATA.chow_10d_hfd_4d_basal ; Raw_data2 = D.RAW_BW_DATA.chow_10d_hfd_4d_end;
elseif kk == 3
Raw_data1 = D.RAW_BW_DATA.chow_8d_hfd_6d_basal ; Raw_data2 = D.RAW_BW_DATA.chow_8d_hfd_6d_end;
elseif kk == 4
Raw_data1 =D.RAW_BW_DATA.hfd_2w_basal ; Raw_data2 = D.RAW_BW_DATA.hfd_2w_end;
end

s = scatter( data.time(1)*ones(length( Raw_data1),1), Raw_data1,'filled',...
'MarkerFaceColor',[ 0 0 0],'SizeData',40);
alpha(s,.45)
s = scatter( data.time(2)*ones(length( Raw_data2),1), Raw_data2,'filled',...
'MarkerFaceColor',[ 0 0 0],'SizeData',40);
alpha(s,.45)
eval(strcat('plot( time_plot, AllSimulations_train.',Experiments_to_plot{kk}, ".variablevalues(:,id_BW),'color', TrainColor - 0.2 , 'linewidth', 2);"))

errorbar(data.time,data.mean, data.SEM,'k*' ,'linewidth',2) % 2wHFD

set(gca,'FontSize',label_size_ax, 'FontWeight','bold')
xlabel('Time (days)','FontSize', label_size )
ylabel('Body weight (g)','FontSize', label_size )
box off
xlim([0 14]);
xticks([0 ,14])
xticklabels({'0','14'})
ylim([22 34]);
ylim([20 34]);

end
%%
%%% figure 3e-f

Experiments_to_plot = {'chow8w','hfd8w'};
for kk = 1:length(Experiments_to_plot)
subplot(2,4,4+kk) %% a)

values=[];
if kk == 1
data.time=Dpred.BW_chow.time; data.mean=Dpred.BW_chow.Mean; data.SEM=Dpred.BW_chow.SEM ;%Control
Expand All @@ -244,17 +265,41 @@

eval(strcat('values(:,1)=Max_',Experiments_to_plot{kk},'_BW;'));
eval(strcat('values(:,2)=Min_',Experiments_to_plot{kk},'_BW;'));


f=fill(time_fill, [ values(1:15,2)' fliplr( values(1:15,1)' ) ], TrainColor, 'EdgeColor', 'none' );
set(f,'facealpha',.5)
hold on
errorbar(data.time(1:3),data.mean(1:3), data.SEM(1:3),'k*' ,'linewidth',1.5) % 2wHFD

fill_time_8w = [ 14:1:56 fliplr(14:1:56) ];

f=fill(fill_time_8w, [ values(15:end,2)' fliplr( values(15:end,1)' ) ], predColor, 'EdgeColor', 'none' );
set(f,'facealpha',.5)
hold on
errorbar(data.time(4:end),data.mean(4:end), data.SEM(4:end),'bO' ,'linewidth',1.5) % 2wHFD

if kk == 1
for kkk = 1:size(D.RAW_BW_DATA.chow_8w_bw,1)
col = [ 0 0 0];
if kkk > 3
col = [ 0 0 1];
end
s = scatter( data.time(kkk)*ones(length( D.RAW_BW_DATA.chow_8w_bw(kkk,:)),1), D.RAW_BW_DATA.chow_8w_bw(kkk,:),'filled',...
'MarkerFaceColor',col,'SizeData',40);
alpha(s,.45)
end
elseif kk == 2

for kkk = 1:size(D.RAW_BW_DATA.hfd_8w_bw,1)
col = [ 0 0 0];
if kkk > 3
col = [ 0 0 1];
end
s = scatter( data.time(kkk)*ones(length( D.RAW_BW_DATA.hfd_8w_bw(kkk,:)),1), D.RAW_BW_DATA.hfd_8w_bw(kkk,:),'filled',...
'MarkerFaceColor',col,'SizeData',40);
alpha(s,.45)
end
end

errorbar(data.time(1:3),data.mean(1:3), data.SEM(1:3),'k*' ,'linewidth',1.5) % 2wHFD
errorbar(data.time(4:end),data.mean(4:end), data.SEM(4:end),'b*' ,'linewidth',1.5) % 2wHFD

% ploting best sim for training and pred
eval(strcat('plot( time_plot, AllSimulations_train.',Experiments_to_plot{kk}, ".variablevalues(:,id_BW),'color', TrainColor - 0.2 , 'linewidth', 2);"))
eval(strcat('plot( time_plot, AllSimulations_pred.',Experiments_to_plot{kk}, ".variablevalues(:,id_BW),'b--', 'linewidth', 2);"))
Expand All @@ -267,18 +312,21 @@
xlim([0 56]);
xticks([0 ,56])
xticklabels({'0','56'})
ylim([22 50]);
ylim([20 50]);

end





%% IPGTT
%% figure 2 - IPGTT - TRAINING AND PREDIKTION
%% figure 4 - IPGTT - TRAINING AND PREDIKTION
%%% figure 2a-d
fill_time_IPGTT = [ t_IPGTT-t_fasting(end) fliplr(t_IPGTT-t_fasting(end)) ]*24*60;
tp_IPGTT=[0 15 30 60 120];
time_best_cost = (t_IPGTT-t_fasting(end) )*24*60;

fig=fig+1;
figure(fig)
set(figure(fig), 'outerposition',[0 0 2560 1440], 'PaperType','a4')
Expand All @@ -302,6 +350,23 @@
set(f,'facealpha',.5)
hold on

if kk == 1
Raw_data = D.RAW_IPGTT_DATA.chow_2w_IPGTT;
elseif kk == 2
Raw_data = D.RAW_IPGTT_DATA.chow_12d_hfd_2d_IPGTT;
elseif kk == 3
Raw_data = D.RAW_IPGTT_DATA.chow_8d_hfd_6d_IPGTT;
elseif kk == 4
Raw_data = D.RAW_IPGTT_DATA.hfd_2w_IPGTT;
end

for kkk = 1:size(Raw_data,1)
col = [ 0 0 0];
s = scatter( data.time(kkk)*ones(length( Raw_data(kkk,:)),1), Raw_data(kkk,:),'filled',...
'MarkerFaceColor',col,'SizeData',40);
alpha(s,.45)
end

errorbar(data.time,data.mean, data.SEM,'k*' ,'linewidth',1.5) % 2wHFD

eval(strcat('plot( time_plot, AllSimulations_train.',IPGTT_sim_name{kk}, ".variablevalues(:,id_Glu),'color', TrainColor - 0.2 , 'linewidth', 2);"))
Expand All @@ -315,29 +380,86 @@

end
%%
%% 6w IPGTT prediction
%%% figure 2e
values=[];
subplot(2,4,5)

hold on
data.time=tp_IPGTT;

data.mean=Dpred.IPGTTgcc.HFDMean; data.SEM=Dpred.IPGTTgcc.HFDSEM;%Control
values(:,1)=Max_hfd6w_Glu; values(:,2)=Min_hfd6w_Glu;

f=fill(fill_time_IPGTT, [ values(:,2)' fliplr( values(:,1)' ) ], predColor ,'EdgeColor', 'none' );
set(f,'facealpha',.5)
hold on

plot( time_plot, AllSimulations_train.IPGTT_hfd6w.variablevalues(:,id_Glu),'color', TrainColor - 0.2 , 'linewidth', 2);
plot( time_plot, AllSimulations_pred.IPGTT_hfd6w.variablevalues(:,id_Glu),'b--' , 'linewidth', 2);

Raw_data = D.RAW_IPGTT_DATA.hfd_6w_IPGTT;
for kkk = 1:size(Raw_data,1)
col = [ 0 0 1];
s = scatter( data.time(kkk)*ones(length( Raw_data(kkk,:)),1), Raw_data(kkk,:),'filled',...
'MarkerFaceColor',col,'SizeData',40);
alpha(s,.45)
end


errorbar(data.time,data.mean, data.SEM,'b*' ,'linewidth',1.5) % 2wHFD
set(gca,'FontSize',label_size_ax, 'FontWeight','bold')
ylabel(' Plasma glucose (mM)','FontSize', label_size ,'FontWeight','bold')
xlabel('Time (min)','FontSize', label_size )
xlim([0 120]);
ylim([5 45]);


%% Basal innsulin
%%% figure 2f
subplot(2,4,6)
values=[];

ins_value_chow2w = AllSimulations_train.IPGTT_chow2w.variablevalues(1,id_Ins) ;
ins_value_chow12dHFD2d = AllSimulations_train.IPGTT_chow12dHFD2d.variablevalues(1,id_Ins) ;
ins_value_chow8dHFD6d = AllSimulations_train.IPGTT_chow8dHFD6d.variablevalues(1,id_Ins) ;
ins_value_hfd2w = AllSimulations_train.IPGTT_hfd2w.variablevalues(1,id_Ins) ;

y = [ D.IPGTTicc.Mean(1,1) ins_value_chow2w ;...
D.IPGTTicc.Mean(1,2) ins_value_chow12dHFD2d ;...
D.IPGTTicc.Mean(1,3) ins_value_chow8dHFD6d ;...
D.IPGTTicc.Mean(1,4) ins_value_hfd2w ];

y = [ NaN ins_value_chow2w ;...
NaN ins_value_chow12dHFD2d ;...
NaN ins_value_chow8dHFD6d ;...
NaN ins_value_hfd2w ];

b=bar(y,'facecolor','flat','HandleVisibility','off');
b(1).CData(1,:) = [ 1 1 1];
b(1).CData(2,:) = [ 1 1 1];
b(1).CData(3,:) = [1 1 1];
b(1).CData(4,:) = [ 1 1 1];

b(2).CData(1,:) = [ 0.7 0.7 0.7];
b(2).CData(2,:) = [ 0.7 0.7 0.7];
b(2).CData(3,:) = [0.7 0.7 0.7];
b(2).CData(4,:) = [ 0.7 0.7 0.7];
hold on

conversion_Factor=0.005;% plasma insulin pmol/L to micro µg/L - 1pmol/L = 0.005µg/L
for kkk = 1:4
if kkk==1
rawdata=D.RAW_IPGTT_DATA.chow_2w_insulin./conversion_Factor;
xt = 0.85;
elseif kkk ==2
rawdata=D.RAW_IPGTT_DATA.chow_12d_hfd_2d_insulin./conversion_Factor ;
xt =1.8500;
elseif kkk ==3
rawdata=D.RAW_IPGTT_DATA.chow_8d_hfd_6d_insulin./conversion_Factor ;
xt =2.8500;
elseif kkk ==4
xt = 3.8500;
rawdata=D.RAW_IPGTT_DATA.hfd_2w_insulin./conversion_Factor ;
end
col = [ 0 0 0];
s = scatter( xt*ones(length( rawdata ),1), rawdata,'filled',...
'MarkerFaceColor',col,'SizeData',40);
alpha(s,.45)
end

errhigh = [Max_chow2w_Ins - ins_value_chow2w,...
Max_chow12dHFD2d_Ins - ins_value_chow12dHFD2d,...
Max_chow8dHFD6d_Ins - ins_value_chow8dHFD6d,...
Expand All @@ -350,9 +472,11 @@


errorbar( 1.15:1:4.15, [ ins_value_chow2w ins_value_chow12dHFD2d ins_value_chow8dHFD6d ins_value_hfd2w],...
errhigh,errlow,'O','color',[ 0.4 0.4 0.4],'linewidth',1.5)
errhigh,errlow,'*','color',[ 0.4 0.4 0.4],'linewidth',1.5)

errorbar( 0.85:1:3.85, [ D.IPGTTicc.Mean(1,1) D.IPGTTicc.Mean(1,2) D.IPGTTicc.Mean(1,3) D.IPGTTicc.Mean(1,4) ], [ D.IPGTTicc.SEM(1,1) D.IPGTTicc.SEM(1,2) D.IPGTTicc.SEM(1,3) D.IPGTTicc.SEM(1,4)],...
'k*', 'linewidth',1.5)
'*', 'color', [0,0,0] , 'linewidth',1.5)

set(gca,'FontSize',label_size_ax, 'FontWeight','bold')
box off
ylabel('Plasma insulin (pM)','FontSize', label_size ,'FontWeight','bold')
Expand All @@ -361,31 +485,5 @@
xtickangle(90)


%% 6w IPGTT prediction
values=[];
subplot(2,4,6)

hold on
data.time=tp_IPGTT;

data.mean=Dpred.IPGTTgcc.HFDMean; data.SEM=Dpred.IPGTTgcc.HFDSEM;%Control
values(:,1)=Max_hfd6w_Glu; values(:,2)=Min_hfd6w_Glu;

f=fill(fill_time_IPGTT, [ values(:,2)' fliplr( values(:,1)' ) ], predColor ,'EdgeColor', 'none' );
set(f,'facealpha',.5)
hold on

plot( time_plot, AllSimulations_train.IPGTT_hfd6w.variablevalues(:,id_Glu),'color', TrainColor - 0.2 , 'linewidth', 2);
plot( time_plot, AllSimulations_pred.IPGTT_hfd6w.variablevalues(:,id_Glu),'b--' , 'linewidth', 2);


errorbar(data.time,data.mean, data.SEM,'bO' ,'linewidth',1.5) % 2wHFD
set(gca,'FontSize',label_size_ax, 'FontWeight','bold')
ylabel(' Plasma glucose (mM)','FontSize', label_size ,'FontWeight','bold')
xlabel('Time (min)','FontSize', label_size )
xlim([0 120]);
ylim([5 45]);




Loading

0 comments on commit 69f73fb

Please sign in to comment.