%%% C:/work/PLPaper2/simul/FinalSumSq.txt

****************************************************************************
******
******  05 July 2005 -- Generate production-quality figures for PLPaper2.
******
******  See C:/work/PLPaper2/simul/NewSpdCmp.txt for simulations with a
******  range of parameter sets. 
******
***************************************************************************
****
****  DEFAULT PARAMETERS for PLM_Hebb2, ver 2.1
****  See C:\work\PLModel1\newsimul\Hebb2m_fits.txt for details.
****
****  For the PLExp1 data set (32 blocks, with feedbk). Based on e32_oN17a.
****    out_gain = 5
****    out_noise = .170    % was .195 in PLPaper1
****    fdbk_wgt = 1.80
****    fdbk_mode = 1 = 'intervene only on error'
****    fdbk_fract = 1 = 'scan all trials'
****    learn_rate = .0016  % was .0015 in PLPaper1 and in ver 2.0
****    criterion = 2.20    % same as in PLPaper1
****    W_init = .17        % same as in PLPaper1 and in ver 2.0
****
****  For the PLExp2 data set (36 blocks, without fdbk). Based on n36_par2.
****  Same parameters as above, except
****    fdbk_fract = 0    (or equivalently, fdbk_mode = 0)
****    out_noise = .156    % was .158 in ver 2.0
****    criterion = 0.95    % was 0.87 in ver 2.0
****
***************************************************************************


>> clear all; cd(fullfile(work_pathstr,'PLPaper2','simul'));
>> global PLM_CACHE;load(fullfile(PLM_pathstr,'PLM_CACHE.mat'));clear Idef;
>> load(fullfile(work_pathstr,'PLPaper2','simul','all_zP_data.mat'));clear mdprfk32 mdprnf36 mZfk32 mZnf36;
>> global HUMAN_ZCONGR ; HUMAN_ZCONGR = eZfk32'; whos

  Name               Size                   Bytes  Class
-------------------------------------------------------------------------
  HUMAN_ZCONGR       6x32                    1536  double array (global)
  PLM_CACHE         12x1                 16802736  struct array (global)
  eZfk32            32x6                     1536  double array
  eZnf36            36x6                     1728  double array
  edprfk32          32x3                      768  double array
  edprnf36          36x3                      864  double array
-------------------------------------------------------------------------

>> Mpar_fk=PLM_params       % defaults, ver 3.3, June 30, 2005
Mpar_fk = 
      rep_size: [7 5 1]
       max_act: 0.5000
      rep_gain: 0.8000
      out_gain: 5
     rep_noise: 0.1000
     out_noise: 0.1700
      W_minmax: [-1 1]
        W_init: 0.1700
       W0_seed: [35x1 double]
    learn_rate: 0.0016
    runav_rate: 0.0200
     criterion: 2.2000
      fdbk_wgt: 1.8000
    fdbk_fract: 1
     fdbk_mode: 1
      blk_size: 300

>> Mpar_nf=Mpar_fk;Mpar_nf.fdbk_fract=0;Mpar_nf.fdbk_mode=0;Mpar_nf.out_noise=.156;Mpar_nf.criterion=.95;Mpar_nf    % PLExp2 defaults
Mpar_nf = 
      rep_size: [7 5 1]
       max_act: 0.5000
      rep_gain: 0.8000
      out_gain: 5
     rep_noise: 0.1000
     out_noise: 0.1560
      W_minmax: [-1 1]
        W_init: 0.1700
       W0_seed: [35x1 double]
    learn_rate: 0.0016
    runav_rate: 0.0200
     criterion: 0.9500
      fdbk_wgt: 1.8000   % irrelevant
    fdbk_fract: 0
     fdbk_mode: 0
      blk_size: 300

>> save(fullfile(work_pathstr,'PLPaper2','simul','FinalSumSq.mat'));


%%%%%%  Run PLM_HEBB2 2000 times on PLExp1 data with feedback params Mpar_fk.
%%
>> N_megaruns=10;N_runs=2*100;tic;for k=1:N_megaruns;sst_fk(k)=NewSpdCmp32(Mpar_fk,N_runs);fprintf('\n');end;toc, 
% Elapsed time 57 minutes on IBM 8310

%% Discard irrelevant sumstat fields and aggregate those four: N_runs, zP_congr, zP_incon, dprime, accZP
>> sst32.N_runs=sum(cat(1,sst_fk(:).N_runs));sst32.zP_congr=mean(cat(3,sst_fk(:).zP_congr),3);sst32.zP_incon=mean(cat(3,sst_fk(:).zP_incon),3);sst32.dprime=mean(cat(3,sst_fk(:).dprime),3);sst32.accZP=mean(cat(1,sst_fk(:).accZP));sst_fk=sst32
sst_fk = 
      N_runs: 2000
    zP_congr: [32x3 double]
    zP_incon: [32x3 double]
      dprime: [32x3 double]
       accZP: [0.8639 0.9483 0.9891 1.3102 0.4811 -0.1462]


%%%%%%  Run PLM_HEBB2 2000 times on PLExp2 data with no-feedback params Mpar_nf.
%%
>> tic;for k=1:N_megaruns;sst_nf(k)=NewSpdCmp36(Mpar_nf,N_runs);fprintf('\n');end;toc, 
% Elapsed time 61 minutes on IBM 8310

%% Discard irrelevant sumstat fields and aggregate those four: N_runs, zP_congr, zP_incon, dprime, accZP
>> sst36.N_runs=sum(cat(1,sst_nf(:).N_runs));sst36.zP_congr=mean(cat(3,sst_nf(:).zP_congr),3);sst36.zP_incon=mean(cat(3,sst_nf(:).zP_incon),3);sst36.dprime=mean(cat(3,sst_nf(:).dprime),3);sst36.accZP=mean(cat(1,sst_nf(:).accZP));sst_nf=sst36
sst_nf = 
      N_runs: 2000
    zP_congr: [36x3 double]
    zP_incon: [36x3 double]
      dprime: [36x3 double]
       accZP: [1.2031 1.3057 1.3581 1.0816 0.1885 -0.4802]

>> clear PLM_CACHE k sst32 sst36;save(fullfile(work_pathstr,'PLPaper2','simul','FinalSumSq.mat'));


%%%%%  Plot Mpar_fk results against the empirical data from PLExp1
>> L=sst_fk;plot_PLM_fit({L.zP_congr',L.zP_incon'},eZfk32');
%% Figure saved as PLPaper2\simul\PLExp1fits_e32_oN17a.pdf    (overwriting the old version)

%%%%%  Plot Mpar(8) results against the empirical data from PLExp2
>> R=sst_nf;plot_PLM_fit2({R.zP_congr',R.zP_incon'},eZnf36');
%% Figure saved as PLPaper2\simul\PLExp2fits_n36_par2.pdf     (overwriting the old version)



***********************************************************
***
***  July 5, 2005 -- Make Figures 8 and 9 (model fits) in PLPaper2/PLFeedback.doc
***  This supersedes earlier versions of Fig 8 & 9 per PLExp2\modelfits\Hebb_fits.txt
***

>> clear all; cd(fullfile(work_pathstr,'PLPaper2','fig'));
>> load(fullfile(work_pathstr,'PLPaper2','simul','FinalSumSq.mat'));


%%  No-feedback d' fits -- top panel of Figure 8 in PLFeedback.doc
 subplot(2,1,1);e=edprnf36;m=sst_nf.dprime;t={[1],[2:9]',[10:17]',[18:25]',[26:33]',[34:36]'};T=[1:36]';
 h=plot(T,e(:,1),'b^',T,e(:,2),'ro',T,e(:,3),'bs',T,m(:,1),'k^',T,m(:,2),'ko',T,m(:,3),'ks',t{2},m(t{2},:),'k-',t{3},m(t{3},:),'k-',t{4},m(t{4},:),'k-',t{5},m(t{5},:),'k-',t{6},m(t{6},:),'k-');
 set(h,'LineWidth',1);set(h(4:6),'MarkerSize',4);   % keep open symbols
 axis([0 37 0 3]);set(gca,'xtick',[0:4:36]);xlabel('Block number  (300 trials/block, 4 blocks/day)');ylabel('d''');title('No feedback  (2000 runs)');

%%  New feedback d' fits in on-error-only mode -- bottom panel of Figure 8 in PLFeedback.doc
 subplot(2,1,2);e=edprfk32;m=sst_fk.dprime;t={[1],[2:9]',[10:17]',[18:25]',[26:31]',[32]};T=[1:32]';
 h=plot(T,e(:,1),'b^',T,e(:,2),'ro',T,e(:,3),'bs',T,m(:,1),'k^',T,m(:,2),'ko',T,m(:,3),'ks',t{2},m(t{2},:),'k-',t{3},m(t{3},:),'k-',t{4},m(t{4},:),'k-',t{5},m(t{5},:),'k-');
 set(h,'LineWidth',1);set(h(4:6),'MarkerSize',4);   % keep open symbols
 axis([0 37 0 3]);set(gca,'xtick',[0:4:36]);xlabel('Block number  (300 trials/block, 4 blocks/day)');ylabel('d''');title('Feedback on error  (2000 runs)');

% Set page size=[1 1 6.5 8 ]; % [Left Bottom Width Height]
% Saved as PLPaper2\fig\Figure08.fig and *.eps
% See PLPaper2\simul\PLExp?fits_*.pdf for the corresponding regression figures

%%  No-feedback incongruent fits -- top panel of Figure 9 in PLFeedback.doc
 subplot(2,1,1);e=eZnf36(:,[2 4 6]);m=sst_nf.zP_incon;t={[1],[2:9]',[10:17]',[18:25]',[26:33]',[34:36]'};T=[1:36]';
 h=plot(T,e(:,1),'b^',T,e(:,2),'ro',T,e(:,3),'bs',t{2},m(t{2},:),'k-',T,m(:,1),'k^',T,m(:,2),'ko',T,m(:,3),'ks',t{3},m(t{3},:),'k-',t{4},m(t{4},:),'k-',t{5},m(t{5},:),'k-',t{6},m(t{6},:),'k-');
 set(h,'LineWidth',1);set(h(7:9),'MarkerSize',4);   % keep open symbols
 axis([0 37 -1 2]);set(gca,'xtick',[0:4:36]);xlabel('Block number  (300 trials/block, 4 blocks/day)');ylabel('z-probability correct');
 title('Incongruent stimuli');axis([0 37 -1 2]);h=refline(0,0);set(h,'Color','k');
 set(gca,'box','on','xtick',[0:4:36],'ytick',[-1:.5:2]);

%%  No-feedback congruent fits -- bottom panel of Figure 9 in PLFeedback.doc
 subplot(2,1,2);e=eZnf36(:,[1 3 5]);m=sst_nf.zP_congr;t={[1],[2:9]',[10:17]',[18:25]',[26:33]',[34:36]'};T=[1:36]';
 h=plot(T,e(:,1),'b^',T,e(:,2),'ro',T,e(:,3),'bs',t{2},m(t{2},:),'k-',T,m(:,1),'k^',T,m(:,2),'ko',T,m(:,3),'ks',t{3},m(t{3},:),'k-',t{4},m(t{4},:),'k-',t{5},m(t{5},:),'k-',t{6},m(t{6},:),'k-');
 set(h,'LineWidth',1);set(h(7:9),'MarkerSize',4);   % keep open symbols
 axis([0 37 -1 2]);set(gca,'xtick',[0:4:36]);xlabel('Block number  (300 trials/block, 4 blocks/day)');ylabel('z-probability correct');
 title('Congruent stimuli');axis([0 37 -1 2]);h=refline(0,0);set(h,'Color','k');
 set(gca,'box','on','xtick',[0:4:36],'ytick',[-1:.5:2]);legend('Data, c=0.245','Data, c=0.160','Data, c=0.106','Model fits',3);

%% Now the fun part --- the inset
 e=[mean(eZnf36(:,[2 4 6])) ; mean(eZnf36(:,[1 3 5]))]
e = 0.9630    0.2839   -0.4317
    1.2077    1.2998    1.3902
 m=[mean(sst_nf.zP_incon) ; mean(sst_nf.zP_congr)]
m = 1.0816    0.1885   -0.4802
    1.2031    1.3057    1.3581

 ah=axes('position',[0.67 0.12 0.22 0.14]);t=[.75 2.25]';T=[1 2];
 h=plot(t,e(:,1),'b^',t,e(:,2),'ro',t,e(:,3),'bs',T,m(:,1),'k^-',T,m(:,2),'ko-',T,m(:,3),'ks-');
 set(h,'LineWidth',1);set(h(4:6),'MarkerSize',4);   % keep open symbols
 axis([.25 2.75 -.75 1.75]);set(gca,'xtick',[1 2],'xticklabel','Incon|Congr','XAxisLocation','top','ytick',[-.5:.5:1.5],'yticklabel','|0||1|');

% Set page size=[1 1 6.5 8 ]; % [Left Bottom Width Height]
% Saved as PLPaper2\fig\Figure09.fig and *.eps


%%%% See PLPaper2\simul\Whist.txt and PLPaper2\simul.Whist.mat for Figure 10


***********************************************************
***
***  July 6, 2005 -- R2 values for the new fits
***

 c=corrcoef(edprfk32(:),sst_fk.dprime(:));c(2)^2
ans =  0.8607      % R2 with feedback. Was .882 in the Psych Review paper

 err=edprfk32-sst_fk.dprime;plot(err,'.-');sqrt(mean(err(:).^2))
ans = 0.2235       % rmse with feedback. Was .209 in the Psych Review paper

 sum(err(:).^2)
ans = 4.7959       % SumSqErr with feedback


 c=corrcoef(edprnf36(:),sst_nf.dprime(:));c(2)^2
ans = 0.9150       % R2 without feedback

 err=edprnf36-sst_nf.dprime;plot(err,'.-');sqrt(mean(err(:).^2))
ans = 0.1821       % rmse without feedback

 sum(err(:).^2)
ans = 3.5815       % SumSqErr without feedback

