2-level PCT model of E. coli

[From Bill Powers (941117.0540 MST)]

Bruce Abbott (various) --

Here is a two-level control system model that works by varying the time
delay before a tumble. The first level is the standard control-system
model Rick and I have been using, in which the delay is gain*error, and
the error is a reference level minus the sensed rate of change of
concentration dNut. The reference level is fixed at 10.

The gain is varied by the output of a second level of control which tries
to keep the sensed concentration at a reference concentration. The error
is multiplied by 0.001 to become the gain of the lower level system (in
other words, control by the second level is through changing the gain of
the first-level system instead of the reference signal as in normal
HPCT).

Changing the fixed reference signal for the first level system primarily
affects the spread in the angles between which the tumbling interval will
be lengthened. If the reference signal RefdNut is set to zero, any
positive dNut will result in lengthening the delay. When RefdNut is set
to 10, dNut must exceed +10 for the delay to be lengthened. Note that
there is a constant 2 added to the timer on each iteration. This creates
a delay of 10 iterations when the error is zero (the trigger level is 20
units). As error goes positive and negative, the delay increases or
decreases around this zero-error level of delay.

The first-level gain determines how sensitive the delay is to a
difference between dNut and RefdNut. The sign of the first-level gain
determines whether E. coli will go up or down the gradient.

If the second-level error RefCon - NewCon is greater than zero, NewCon is
too small, and gain is positive. A positive gain causes E. coli to swim
down the gradient (away from the target) because when dNut is less than
RefdNut, "timer" is increased more rapidly on each iteration and the
delay to the next tumble is shortened.

Note that it is the setting of the second-level reference signal that
determines what level of concentration the system will seek -- that is,
what level of concentration it appears to find "rewarding."
Concentrations both higher and lower than that level are observed to be
less "reinforcing" or even "punishing" (the system will avoid
concentrations higher or lower than the reference level). The reinforcing
or punishing quality, of course, does not lie in the sensed
concentration, because the system will seek any concentration specified
by the reference signal RefCon inside the control system, and in doing so
will appear to avoid other levels of concentration as if they were
aversive. When the reference signal RefCon is decreased, what had seemed
to be a rewarding level of concentration suddenly becomes aversive -- the
system flees outward from the target, until it finds the specified lower
level of concentration.

Since the variable being sensed changes only with radial distance from
the target, the controlled consequence of behavior is only this radial
distance. Random movements around a circle of constant radius are not
resisted; they have no effect on sensed concentration. As you change
RefCon with the + and - keys, you see the radius of this circle changing.

From the numbers displayed, you can see that the sensed concentation

NewCon is always kept near the reference concentration RefCon.

Best,

Bill P.

{********************************************************************
* E. COLI "REINFORCEMENT" SIMULATION 4awp *
* *
* *
* Modification of Bruce Abbott's Ecoli4a to use continuous control *
* and to seek a specific level of concentration by varying the gain *
* of the loop that controls rate of change of concentration. *
* *
* Change ref level of concentration with + and - keys in ratio *
* steps of 1.1. *
* *
* WTP, 941117 *
* *
***********************************************************************}

program Ecol4awp;
uses
  CRT, Graph, GrUtils;

const
  TWOPI = PI * 2;
  ENDSESSION = 50000;
var
  MaxX, MaxY: integer;
  NutrX, NutrY, X, Y: integer;
  dNut,NutMag,OldCon,NewCon,RefCon: real;
  EcoliX, EcoliY, timer,timelim, RefdNut, gain, error,
  Speed, Angle: real;
  Ch: char;
  Clock: longint;
  t: integer;

procedure InitScreen;
begin
  ClrScr;
  InitGraphics;
  MaxX := GetMaxX; MaxY := GetMaxY;
  Rectangle(0, 0, MaxX, MaxY);
  OutTextXY(MaxX div 2 - 220, Y+5,
    'E. COLI SIMULATION # 4AWP: 2-LEVEL CONTROL MODEL');
  OutTextXY(MaxX - 200, 50, 'C ref level');
  OutTextXY(MaxX - 200, 60, 'C perceived ');
  OutTextXY(20, MaxY-50, 'Press ESC or q to Quit...');
  OutTextXY(20, MaxY-35, 'Press + or - to change ref concentration...');
  OutTextXY(20, MaxY-20, 'Press space to pause or continue...');
end;

procedure ShowReal(x,y: integer; v: real);
var s: string;
begin
str(v:6:1, s);
setfillstyle(0,0);
bar (x,y,x+textwidth(s),y+textheight(s));
outtextxy(x,y,s);
end;

procedure InitSim;
begin
  Randomize;
  Speed := 1.0;
  NutrX := MaxX div 2;
  NutrY := MaxY div 2;
  EcoliY := 10 + random(100);
  EcoliY := 10 + random(100);
  Angle := TwoPi * Random;
  EcoliX := EcoliX + Speed * cos(Angle);
  EcoliY := EcoliY + Speed * sin(Angle);
  X := Round(maxx div 2 + EcoliX);
  Y := Round(maxy div 2 - EcoliY);
  Rectangle(NutrX-2, NutrY-2, NutrX+2, NutrY+2);
  NutMag := 100.0; { max concentration }
  Clock := 0;
  timelim := 20.0;
  timer := 0;
  RefdNut := 10;
  RefCon := 5000;
  gain := -0.01;
end;

{ADJUSTED TO GIVE LARGER VALUES OF NUTCONCEN}
function NutConcen(X, Y: real): real;
{ Nutient concentration at point X, Y: environment function }
var
  Dist: real;
begin
  Dist := Sqr(X) + Sqr(Y);
  NutConcen := 3E5*NutMag / (100.0 + 0.001*(Sqr(Dist)));
end;

procedure StepEColi;
var
  NewNut: real;
begin
  EcoliX := EcoliX + Speed * cos(Angle); {MOVE E COLI}
  EcoliY := EcoliY + Speed * sin(Angle);
  X := Round(EcoliX);
  Y := Round(EcoliY);
  PutPixel(maxx div 2 + X, maxy div 2 - Y, white); {PLOT POSITION}
  NewCon := NutConcen(EcoliX,EcoliY); {GET CONCENTRATION}
  dNut := NewCon - OldCon; {CALCULATE CHANGE}
  OldCon := NewCon;
  {GAIN OF LEVEL 1 SYSTEM = 0.001 * ERROR IN SECOND LEVEL SYSTEM}
  gain := 0.001*(RefCon - NewCon);
  timer := timer + 2 + gain * (RefdNut - dNut); {STEP TOWARD TUMBLE]
  if timer < 0.0 then timer := 0.0;
  if timer > timelim then
   begin
    Angle := TwoPi * Random; { TUMBLE}
    timer := 0.0; { RESET TIMER}
   end;
  showreal(MaxX - 100, 50,RefCon);
  showreal(MaxX - 100, 60,NewCon);

end;

begin
  InitScreen;
  InitSim;
  randomize;
  repeat
    inc(Clock);
    ch := chr(0);
    StepEcoli;
    if Keypressed then Ch := ReadKey;
    if ch = '+' then if RefCon < 50000.0 then RefCon := RefCon*1.1;
    if ch = '-' then if Refcon > 500.0 then RefCon := RefCon/1.1;
    if ch = ' ' then ch := readkey;
  until (Ch = #27) or (Clock >= ENDSESSION) or (ch = 'q');
  RestoreCRTMode;
  CloseGraph;
end.