%
% hopstat.m
%
% train and run a hopfield dna ocr net
%

function hopstat(nrun)

clf

P = -ones(25,4);
P([1:5 8 13 18 23],1) = 1; 	  	% T
P([1:5 6 10 11:15 16 20 21 25],2) = 1;	% A
P([1:5 6 11 16 21:25],3) = 1;		% C 
P([1:5 6 11 14 15 16 20:25],4) = 1;	% G 

W = P*P';

hit = zeros(1,9);

s = sign2(rand(25,1));

if nrun==1,
   subplot(1,2,1)
   led(s)
end

for n=1:nrun

    err = 1;

    while err>0

        ns = sign2(W*s);
    
        err = max(abs(s-ns));

        s = ns;

    end

    hit
    led(s)

    [tf,loc] = ismember(s',P','rows');

    if tf,
       hit(loc) = hit(loc)+1;
    else
       [tf,loc] = ismember(-s',P','rows');
       if tf,
          hit(4+loc) = hit(4+loc)+1;
       else
          hit(9) = hit(9)+1;
       end
    end
    hit
    pause

    s = sign2(randn(25,1));
       
end

hit = hit/nrun;

if nrun==1,
   subplot(1,2,2)
   led(ns)
   return
end
sum(hit)
bar(hit)
lab = {'T', 'A', 'C', 'G', '-T', '-A', '-C', '-G', '?'};
set(gca,'xticklabel',lab)
xlabel('pattern','fontsize',16)
ylabel('probability','fontsize',16)
title(['Hopfield OCR response to ' num2str(nrun) ' stimuli'],'fontsize',16)

function led(state)
for j=1:25,
    y = fix((j-.1)/5);
    x = j-5*y-1;
    if state(j) > 0
       c = [1 0 0];
    else
       c = [1 1 1];
    end
    patch([x x+1 x+1 x],[-y -y -y-1 -y-1],c)
end
axis([0 5 -5 0])
axis off
axis equal
drawnow


function val = sign2(x)
tmp = sign(x);
val = tmp + abs(tmp) - 1;