forked from mpc-ucl/emfit
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathbatchParameterPlots.m
93 lines (76 loc) · 2.37 KB
/
batchParameterPlots.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
function batchParameterPlots(Data,models,resultsDir,bestmodel)
nModls = length(models);
Nsj = length(Data);
nfig=0;
eval(['mkdir ' resultsDir '/figs']);
%--------------------------------------------------------------------
% plot parameters of best model
%--------------------------------------------------------------------
figure(2);clf;
for mdl=1:nModls
R.(models(mdl).name) = load(sprintf('%s/%s',resultsDir,models(mdl).name));
end
for j=1:3
if j==1
% posterior parameters
par = R.(models(bestmodel).name).E;
yl = models(bestmodel).parnames_untr;
yl{1} = {'EM-MAP',yl{1}};
elseif j==2
% transformed posterior parameters
R = modelConvertParams(R,models);
par = R.(models(bestmodel).name).parConvert;
yl = models(bestmodel).parnames;
yl{1} = {'EM-MAP',yl{1}};
elseif j==3
% lightly regularized ML parameters - closest to the data
par = R.(models(bestmodel).name).stats.EMAP0;
npar= models(bestmodel).npar;
yl = models(bestmodel).parnames_untr;
yl{1} = {'MAP',yl{1}};
end
npar = size(par,1);
m = mean(par,2);
s = std(par,[],2)/sqrt(Nsj);
for k=1:npar
subplot(3,npar,k+(j-1)*npar)
mybar(m(k),.7);
hon
plot(linspace(0.7,1.3,Nsj),par(k,:),'k.');
mydeb(0,m(k),s(k));
hof
set(gca,'xtick',[]);
ylabel(yl{k});
end
end
myfig(gcf,[resultsDir filesep 'figs/Parameters']);
%--------------------------------------------------------------------
% try plotting parameters against true values if ran on generated data
%--------------------------------------------------------------------
if isfield(Data,'trueModel');
try
parTrue = [Data.trueParam];
parFitEM = R.(Data(1).trueModel).E;
parFitML = R.(Data(1).trueModel).stats.EML;
truemodelid = find(strcmp({models.name},Data(1).trueModel));
npar = size(parTrue,1);
figure(3);clf;
for k=1:npar
subplot(1,npar,k)
plot(parTrue(k,:),parTrue(k,:),'r-','linewidth',2);
hon
h(1) = plot(parTrue(k,:),parFitEM(k,:),'b+','markersize',20);
h(2) = plot(parTrue(k,:),parFitML(k,:),'k.','markersize',20);
hof
xlabel('True Param Value');
if k==1; ylabel('inferred param value');end
if k==1; legend(h,{'EM-MAP','ML'},'location','northwest');end
title(models(truemodelid).parnames_untr(k));
set(gca,'fontsize',14);
mytightaxes
end
myfig(gcf,[resultsDir filesep 'figs/ParameterRecovery']);
catch
fprintf('Plotting true against inferred parameters failed.\n');
end
end