function  F = PLE_2exp(params,xdata)

%PLE_2EXP -- Perceptual learning: two parallel exponential curves.
%
%  F = PLE_exponent(params,xdata)
%  F = PLE_exponent(params)
%
%  This function is designed to be minimized via PLE_CURVEFIT or LSQCURVEFIT.
%  It implements two parallel exponential curves, exploring the
%  possibility that each context learns independently of the other.
%  Two nested models are explored: with 3 and 5 free params, resp.
%  The length of the PARAMS argument determines which formula is used:
%  Model 1:  F(R) = F_as + (F_in-F_as)*exp(-alpha*(R-1))
%    where R grows independently for contexts A and B.
%  Model 2:  F(T,R) = F_as + (F_in-F_as)*exp(-alpha*(T-1)) + ...
%                   + G_as        -G_as *exp(-beta*(R-1))
%    where T=[1:32] and R grows independently for contexts A and B.
%
%  PARAMS is a vector of length either 3 or 6 specifying parameters:
%   1. F_as -- asymptotic value of the general curve (when T=Inf)
%   2. F_in -- initial value of the general curve (when T=1)
%   3. alpha-- learning rate of the context-general mechanisms
%   4. G_as -- asymptotic value the specific curve (when R=Inf)
%   5. beta -- learning rate of the context-specific mechanisms
%
%  XDATA depends on the model. For model 1, it is a 32x1 column vector
%  specifying R in Eq.1 above. Default [1,1:8,2:9,9:16,10:15,17]'.
%  For model 2, XDATA is a 32x2 matrix specifying T and R in Eq.2.
%  Default [1:32 ; 1,1:8,2:9,9:16,10:15,17]'.
%
%  Two parallel learning curves, Model 1:
%   edpr=mean(HUMAN_ZCONGR)'.*2;blk1=[1,10:17,26:31]';blk2=[2:9,18:25,32]';
%   x0=[2 1 .2];T2=[1,1:8,2:9,9:16,10:15,17]';
%   opt_par=lsqcurvefit('PLE_2exp',x0,T2,edpr),mdpr=PLE_2exp(opt_par);
%   plot(1:32,mdpr,'k.-',blk1,edpr(blk1),'ro',blk2,edpr(blk2),'bs');
%
%  Null hypothesis: single overarching curve, Model 1:
%   T1=[1:32]';opt_par1=lsqcurvefit('PLE_2exp',x0,T2,edpr)
%   plot(1:32,PLE_2exp(opt_par1,T1),'k.-',1:32,edpr,'ro');
%
%  See also PLE_EXPONENT, PLE_CURVEFIT, LSQCURVEFIT.

%  Original coding by Alexander Petrov, UC Irvine.
%  $Revision: 1.0 $  $Date: 2002/12/12 14:26 $

%-- Decipher arguments
if (length(params)==3)         % Model 1
   if (nargin < 2)
      F = model1(params,[1,1:8,2:9,9:16,10:15,17]') ;
   else
      F = model1(params,xdata) ;
   end
elseif  (length(params)==5)    % Model 2
   if (nargin < 2)
      F = model2(params,[1:32]',[1,1:8,2:9,9:16,10:15,17]') ;
   else
      F = model2(params,xdata(:,1),xdata(:,2)) ;
   end
else
   error('PARAMS must have either 3 (Model 1) or 5 (Model 2) elements.') ;
end

%--  Return F
%--- End of the main function


%%%%%  Subfunctions  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function  F = model1(params,R)

F_as = params(1) ;                 % asymptotic level
delta_F = params(2) - F_as ;       % improvement = F_initial - F_asympt
minus_alpha = -params(3) ;         % learning rate

F = F_as + delta_F.*exp(minus_alpha.*(R-1)) ;   % Equation 1
%--  Return F
%--- End of subfunction MODEL1  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function  F = model2(params,T,R)

F_as = params(1) ;                 % general asymptotic level
delta_F = params(2) - F_as ;       % improvement = F_initial - F_asympt
minus_alpha = -params(3) ;         % general learning rate
G_as = params(4) ;                 % specific asymptotic level
minus_beta = -params(5) ;          % specific learning rate

F = F_as + delta_F.*exp(minus_alpha.*(T-1)) + ...
    G_as -    G_as.*exp(minus_beta.*(R-1))  ;        % Equation 2
%--  Return F
%--- End of subfunction MODEL2

%%%%% End of PLE_2EXP.M
