testrtrl.m

two Matlab functions for initializing and training a recurrent neural network

choice = input(...
    [' 1. 2-bit xor, delay is 2, 1 cycles\n' ...
      ' 2. 1-bit identity, 4 step delay, 1 cycle, random input\n' ...
      ' 3. bit 1 followed indefinitely by bit 2, 1 cycle, random input\n' ...
      ' 4. second order system, given state variables as input\n' ...
      ' 5. second order system, given history of one state variable as input\n' ...
      ' 6. third order system, suspension system, given one state variable' ...
      ' and action as input \n' ...
      '    Enter choice: ']);
if choice == 1

  figure(1);
  ni = 2; no = 1; nu = 4; cycles=1; delay = 2;
  net = rtrlinit(ni,no,nu,'sigmoid',0.1);
  state = zeros(1,delay*2);
  numsamples = 100000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  for isample = 1:numsamples
    targ = xor(state(1),state(2))*0.8+0.1;
    x = rand(1,2)>0.5;
    sample = [x targ];
    state = [state(3:end) x]; 
    net = rtrl(net,cycles,sample,1);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

elseif choice == 2

  figure(1);
  ni = 1; no = 1; nu = 4; cycles=1; delay = 4; alpha = 1;
  net = rtrlinit(ni,no,nu,'sigmoid',0.1);
  state = zeros(1,delay);
  numsamples = 40000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  for isample = 1:numsamples
    x = rand(1,1)>0.5;
    targ = state(1) *0.8+0.1;
    sample = [x targ];
    state = [state(2:end) x];
    net = rtrl(net,cycles,sample,alpha);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

elseif choice == 3

  figure(1);
  ni = 2; no = 1; nu = 2; cycles=1; delay = 0; alpha = 2;
  net = rtrlinit(ni,no,nu,'sigmoid',0.1);
  bit1appeared = 0;
  numsamples = 50000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  for isample = 1:numsamples
    x = rand(1,ni)>0.8;
    bit1appeared = bit1appeared | x(1);
    targ = (bit1appeared & x(2))*0.8+0.1;
    if targ > 0.5
      bit1appeared = 0;
    end
    sample = [x targ];
    net = rtrl(net,cycles,sample,alpha);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

elseif choice == 4

  figure(1);
  ni = 2; no = 1; nu = 2; cycles=1; delay = 0; alpha = 0.1;
  net = rtrlinit(ni,no,nu,'tanh',0.1);
  state = zeros(1,(delay+1)*ni);
  numsamples = 10000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  %A = [0 -2.5; 4 -2];  B = [0 -1.25 ; -0.5 -0.5]; C = [0 2]; D = [0 1];
  x = [0;0.1];
  for isample = 1:numsamples
    %u = (rand(2,1)-0.5)*0.1
    x = x + 0.05 * ([0 1; -2 0] * x - [0 0.1] * x);
    %y = C * x + D * u
    targ = x(1);
    sample = [state(1:2) targ] * 10;
    state = [state(3:end) x'];
    
    net = rtrl(net,cycles,sample,alpha);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

elseif choice == 5

  figure(1);
  ni = 2; no = 1; nu = 3; cycles=1; delay = 0; alpha = 0.1;
  net = rtrlinit(ni,no,nu,'sigmoid',0.1);
  state = zeros(1,2);
  numsamples = 20000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  %A = [0 -2.5; 4 -2];  B = [0 -1.25 ; -0.5 -0.5]; C = [0 2]; D = [0 1];
  x = [0;0.1];
  for isample = 1:numsamples
    %u = (rand(2,1)-0.5)*0.1
    x = x + 0.05 * ([0 1; -2 0] * x - [0 0.1] * x);
    %y = C * x + D * u
    targ = x(1);
    sample = ([state(1:2) targ]+0.1) * 5;
    state = [state(2:end) x(1)];
    
    net = rtrl(net,cycles,sample,alpha);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

elseif choice == 6
  
  %http://www.engin.umich.edu/group/ctm  
  figure(1);
  ni = 2; no = 1; nu = 5; cycles=1; delay = 0; alpha = 0.1;
  net = rtrlinit(ni,no,nu,'sigmoid',0.1);
  state = zeros(1,1);
  numsamples = 200000;
  samplesperplot = max([1 numsamples/100]);
  sumerr = 0; errs=[]; outputs = zeros(100*cycles,no+nu+ni);
  A = 1.e3 * [       0    0.0010         0         0
   -0.0066         0   -0.0253   -0.0001
    0.0469         0   -0.0482    0.0010
    1.5625         0   -1.8445         0];
  B =1.e3 * [0         0
    0.0000    0.0066
         0   -0.0469
    0.0000   -1.5625];
  %Calc norm parms
  x = [0;0;0;0];
  n = 1000;
  samples = zeros(n,3);
  for is = 1:10000
    u = randn(1,1)*1;
    dx =  A * x + B * [u;1];
    x = x + 0.01 * dx;
    samples(is,:) = [x(3) u dx(3)];
  end
  means = mean(samples)
  vars = var(samples)
  mins = min(samples)
  maxs = max(samples)
  %Now run it
  x = [0;0;0;0];
  us = [];
  for isample = 1:numsamples
    u = u - 0.1 * (u + randn(1,1));
    dx =  A * x + B * [u;1];
    us = [us dx];
    x = x + 0.01 * dx;
    %y = C * x + D * u
    targ = x(3);
    sample = [([state(1) u]-means(1:2))./vars(1:2)*3  (dx(3)-mins(3))/(maxs(3)-mins(3))];
    
    state = x(3);
    net = rtrl(net,cycles,sample,alpha);
    [errs,sumerr,outputs] = testrtrlplots(isample,samplesperplot,...
	cycles,errs,outputs,sumerr,sample,net);
  end

end