% TP6ex1 : le bandit à deux bras

clear
clf
A=input('Probabilité de gain associé au levier A : ');
B=input('Probabilité de gain associé au levier B : ');
n=input('Nombre d''itérations : ');

% initialisation des variables
NA=0;NB=0;TA=0;TB=0;
G=[];
ESTA=[];
ESTB=[];
%variables supplémentaires pour la stratégie optimale
d=2;
p=round(log(n)/log(d));
PA=[2:2:p];
PB=[1:2:p];
IA=2.^PA;
IB=2.^PB;
% on initialise le processus en choisissant avec une chance sur deux l'un des deux leviers.
z=rand;
if z<0.5
   U=A;NA=1;
else 
   U=B;NB=1;
end

% Debut du programme
for i=1:n
   % Calcul de la sortie à l'instant i
   r=rand;
   if r<(1-U)
      X=0;
   else
      X=1;
   end
   % Estimateurs des moindres carrés à l'instant i
   if (U==A & X==1)
      TA=TA+(1-TA)/NA;
   elseif (U==A & X==0);
      TA=TA-TA/NA;
   else
      TA=TA;
   end
   if (U==B & X==1)
      TB=TB+(1-TB)/NB;
   elseif (U==B & X==0);
      TB=TB-TB/NB;
   else
      TB=TB;
   end
   
   % Choix du levier à l'instant i : mauvaise stratégie (quoique naturelle)
%   if TA>=TB
 %     U=A;NA=NA+1;
  % else
   %   U=B;NB=NB+1;
   %end
   
   % Choix du levier à l'instant i : stratégie optimale
   hA=ismember(i,IA);
   hB=ismember(i,IB);
   if hA==1
      U=A;NA=NA+1;
   elseif hB==1
      U=B;NB=NB+1;
   elseif TA>=TB
      U=A;NA=NA+1;
   else
      U=B;NB=NB+1;
   end
   % Stockage
   G=[G X];
   ESTA=[ESTA TA];
   ESTB=[ESTB TB];
   % fin de la boucle for
end

% Sorties graphiques

subplot (3,1,1)
plot(cumsum(G)./[1:length(G)],'r')
title('Stratégie optimale pour le bandit à deux bras')
ylabel('Gains moyens')
hold on
C=max(A,B);
plot(C*ones(n,1),'b--')
grid

subplot(3,1,2)
plot(ESTA,'r')
ylabel('Estimation de A')
hold on
plot (A*ones(n,1),'b--')
grid

subplot(3,1,3)
plot(ESTB,'r')
ylabel('Estimation de B')
xlabel('Nombre de réalisations')
hold on
plot (B*ones(n,1),'b--')
grid

hold off

% Je suis tombé dans un piège avec les valeurs A=0.6,B=0.5,N=10000