------------------- A.A follows --------------------
[Hans Blom, 950530]
(Bill Powers (950523.2115 MDT)), the "missing post"
I don't understand the Kalman approach, but perhaps you can figure this
out. Assuming that the world-model has converged to match the real
system (a = at, b = bt, c = ct), can you solve the equations
analytically for u?
When the world-model has fully converged (a=at etc), several things must
be noted. First, the prediction by the model is
x [k+1] = ct + at * x [k] + bt * u [k]
whereas the "world" acts according to
xt [k+1] = ct + at * xt [k] + bt * u [k] + ft [k]
where ft stands for the unmodelled dynamics. Second, the control
algorithm that I implemented in my demo simply attempts to compute a
u [k] such that the predicted x is taken to xopt in one step, i.e.
xopt = ct + at * x [k] + bt * u [k]
or
u [k] = (xopt - ct - at * x [k]) / bt
Third, a fully converged model implies that all information from the
world has been fully extracted, and that the observation y plays no role
whatsoever anymore: you don't have to look anymore, because you know. So
My suspicion is that u will turn out to be exactly a constant times x +
disturbance, where disturbance = xt - x.
is incorrect, as well as all speculation based upon it. Note that we are
talking about a fully converged model here, in which "forgetting" plays
no role.
RE: a counterexample to control of perception
My demo includes a noise term (vt) that corrupts the perception
only, NOT what happens in the world (that is modelled by the term
ft). Set vt to a non-zero value, and see what is controlled: the
world (xt), NOT the perception (y).
That is true under the assumption that the variations in the perceptual
signal y are not true representations of xt. Remember that your noise
term is a high-frequency noise; your adaptation process can't keep up
with it, and effectively smooths it out.
Frequency plays no role. In the code that follows I present a demo where
the observation noise term vt contains both a low frequency AND a high
frequency component (ft = 0). After convergence, both are effectively
disregarded so that xt tracks xopt, whereas y does NOT track xopt. Note
that, after convergence, the position of the y-trace may, for long
periods, lay either above or below xopt. That does not matter.
If the same variations in xt
were in fact real (i.e., caused by some other noise contribution), and
the perceptual function were noiseless, your system would still ignore
the fast variations and control only the low-frequency behavior of xt.
Again, frequency plays no role. What would follow xopt is x, the modelled
part of xt. Whether the unmodelled dynamics is high- or low-frequency is
not important.
(Bill Powers (950529.1230 MDT))
What is new in my demo is, in comparison to the "standard"
PCT-controllers:
1) an internal model of the regularities of "the world out there"
can be built, so that for instance a square wave pattern reference
level can be accurately tracked with essentially no error and zero
response delay;
If you're speaking of the model you presented, I don't think this is
true. The optimum fit of world-model to external system would result in
control having exactly the same accuracy and the same delay as direct
control without the world-model. I think you had better check this claim
against your model run with a faster time-scale.
Remember that in a fully converged model the observations are super-
fluous: you do not need to check whether you are doing it right, you KNOW
it. Run my demo with a reference level xopt [*] of as high a frequency
sine or square wave as you like, and verify that for all frequencies
xt [*] will track xopt [*] with no delay.
3) this internal model may become so trustworthy that it is more relied
on than the perceptions; this, as well as 2) above, was my counter-
example where the generalization "control of perception" breaks down
and must be replaced by "control of (by?) internal representations".
You're making an unwarranted assumption here, that the world-model is
more trustworthy than the perceived state of the actual system. I
believe that this assumption, whether in this particular model or in
real life, is mistaken more often then correct.
I say "this internal model MAY become so trustworthy that ..." I do not
imply that we will always be able to build such a trustworthy model, only
that SOMETIMES we will. You seem to agree to that.
When you proposed the demonstration refuting control of perception, you
used a disturbance consisting of high-frequency random noise.
Subsequent speculations are out of order. Frequency plays no role. See
the code of the demo that follows.
There is no way to distinguish the variations in y due to the
disturbance from variations due to xt.
This is where your misunderstanding of Kalman filtering shows up most. In
fact, the Kalman filter IS able to make this distinction. That is the
reason why it has TWO noise sources, ft (the unpredictable part of xt)
and vt (sensor noise), which are treated quite differently. Verify this
in the following demo and check the different behavior for ft = 0 and vt
<> 0 versus ft <> 0 and vt = 0. In particular, note that
The disturbance of y is actually causing the real xt to vary. The result
is that y, but not xt, is following the changes in xopt.
is incorrect.
Nothing of this is new, of course. 1) Tracking a regular
("predictable") pattern with approximately zero delay has
frequently been encountered in human operator response speed
studies.
You will find a demonstration of this effect in the 1960 paper by
myself, Clark, and MacFarland. This can be be explained, as I said,
using a hierarchical model.
Yes, I am sure that different models can "explain" (mimick? is that what
you mean?) the same observations.
... The value of x is always being affected by x - y, even after
convergence is complete and the coefficients are no longer changing.
Run the demo and test whether this is true and under which conditions...
The system noise, being zero on average, is disregarded.
It is only through disregarding the system noise that you can come to
the conclusion that the model-based system is accurately representing
the real system.
No, the model-based system is accurately representing the real system
ONLY IF the system noise is exactly zero. Only in that case do model and
world coincide. On the other hand, the model-based system is THE BEST
POSSIBLE representation of the real system if no properties the system
noise can be modelled anymore, i.e. when the system noise is white noise.
But this is valid only for the high-frequency components of the system
noise. If you were to do a frequency analysis of the behavior of the
model, you would find that the low-frequency components are being
reduced by the loop gain of the system, and xt is actually being made
to vary at these frequencies in a way that x does not vary.
Run the demo below and verify that this is not true.
Again, it is not only the ability to operate blind that I
find interesting. It is the ability to somehow store
knowledge about the world that we can use in our control
tasks IN ADDITION TO our direct perceptions.
Yes, that is an interesting subject. But it is of interest in studies of
human behavior only if we can show that there are real behaviors that
can't be accounted for economically in any simpler way.
Of interest... To whom? We can look at objects in terms of atoms, but we
can also chunk in many different ways. But let's not get stuck in more
philosophy.
And what we can observe is restricted, as you note, by what
we can accurately measure, i.e. in terms of what is already
known to us. This embeds observations within a culture: the
type of observations that we can make will always depend on
whichever instruments we have and on our pre-existing
knowledge of what is observation-worthy.
Well, yeah, I suppose so, but so what?
Well, yeah, what for you is a "so what" is for me of fundamental
importance, because it shows me my limitations -- and those of others.
That which you have learned to see and to be worthy of attention is
somewhat like your personally discovered set of axioms on which you build
your theories and which determine what you (can) achieve in life. Your
private axioms are different from those of others, so others will arrive
at different conclusions and lead different lifes. Others will only be
able to arrive at your conclusions if you teach them your axioms, i.e. if
you convince them to see what you see, to find important what you find
important, to lead your sort of life. And I bet that, except for a few
devout followers, few are willing to do just that. In the above, replace
"set of axioms" by "world-model" to see where my interests lie.
I believe that your Kalman-filter program achieves resistance to
external disturbances because there is an unforseen connection between
the perceptual input xt and the comparator (where xopt comes into play).
This unforseen connection gives the true value of xt a direct effect on
the signal entering the comparator, and thus completes a negative
feedback loop that is largely independent of the adaptation processes.
An "unforeseen" connection? Magic? That reminds me of a saying:
Science is magic to those who do not know its rules.
Magic becomes science for those who learn its rules.
Isn't that basically about building an ever better world-model with ever
less unmodelled dynamics?
Here is the code to which I referred above. For easier experimentation I
have left out those parts that deal with parameter changes and forget-
ting, so that convergence can be reached. I have also moved the definit-
ions of the variables/constants ft, pff, vt, and pvv to the top of the
program, renaming them to fto, pffo, vto, and pvvo. This new demo starts
out with fto=pffo=0 (no unmodelled dynamics) and vt a low frequency plus
a high frequency sinusoidal waveform. Check that both the low as well as
the high frequency components in the observation noise are going to be
eliminated and that xt will start to track xopt, regardless of the value
of y. Check also that after convergence (after run 3000; actually pseudo-
convergence, because the parameter uncertainty is not yet equal to zero
and the parameters have not yet reached their true values) the observation
is not required anymore.
Demo: Pseudo-control: Control of INNER, not OUTER perception.
program kalman_filter_simulation_with_two_sines;
{directions for use of the program:
compile with Turbo Pascal and run
(must have EGAVGA.BGI in current directory!)
press SPACE to view the next batch of runs
press RETURN to end simulation}
uses
crt, graph;
const
fto = 0.0;
pffo = 0.000001;
vto = 0.5;
pvvo = 0.5; {TWO sine waves of amplitude vto}
minX = 40; {start of plots at this x}
var
run: longint; {increments each sample}
batch: integer; {increments each batch}
t: integer; {step counter within batch}
{system variables}
xt, at, bt, ct: real; {system parameter values}
ft: real; {system noise amplitude}
vt: real; {observation noise amplitude}
xopt: real; {optimal x, setpoint}
{model variables}
x, a, b, c: real; {model parameter values}
pxx, pxa, pxb, pxc,
pax, paa, pab, pac,
pbx, pba, pbb, pbc,
pcx, pca, pcb, pcc: real; {and their covariance matrix}
pff: real; {assumed system noise variance}
pvv: real; {assumed observ. noise variance}
u: real; {action}
y: real; {observation}
r: real; {residual = prediction error}
pyy: real; {expected variance of residual}
qaa, qbb, qcc: real; {expected variation on change}
{screen variables}
MaxX, MaxY: integer; {screen dimensions}
xmin, xmax,
amin, amax,
bmin, bmax,
cmin, cmax,
umin, umax: real; {display scaling factors}
{procedures for plotting}
procedure var_scales;
{define minima and maxima for plotting; change if you like}
begin
xmin := 8.0; xmax := 12.0;
amin := 0.8; amax := 1.0;
bmin := -0.2; bmax := 0.2;
cmin := 0.0; cmax := 2.0;
umin := -6.0; umax := 6.0
end;
procedure init_graphics; {initialize graphics}
var
grDriver,
grMode,
ErrCode,
i, j: integer;
s: string [50];
begin
grDriver := Detect;
InitGraph (grDriver, grMode, '');
ErrCode := GraphResult;
if ErrCode <> grOk then
begin
WriteLn('Graphics error:', GraphErrorMsg (ErrCode));
halt
end;
setcolor (white);
{hope this works on other than EGA/VGA displays!}
maxX := GetMaxX; maxY := GetMaxY - 20;
line (minX, 0, minX, maxY);
line (minX, maxY, maxX, maxY);
line (maxX, maxY, maxX, 0);
for i := 0 to 4 do
begin
line (minX, i * maxY div 5, maxX, i * maxY div 5);
moveto (0, i * maxY div 5 + maxY div 10);
case i of
0: outtext (' x');
1: outtext (' a');
2: outtext (' b');
3: outtext (' c');
4: outtext (' u');
end
end;
for i := 0 to 4 do
begin
moveto (0, i * maxY div 5 + 5);
case i of
0: begin str (xmax:4:1, s); outtext (s); end;
1: begin str (amax:4:1, s); outtext (s); end;
2: begin str (bmax:4:1, s); outtext (s); end;
3: begin str (cmax:4:1, s); outtext (s); end;
4: begin str (umax:4:1, s); outtext (s); end;
end
end;
for i := 1 to 5 do
begin
moveto (0, i * maxY div 5 - 10);
case i of
1: begin str (xmin:4:1, s); outtext (s); end;
2: begin str (amin:4:1, s); outtext (s); end;
3: begin str (bmin:4:1, s); outtext (s); end;
4: begin str (cmin:4:1, s); outtext (s); end;
5: begin str (umin:4:1, s); outtext (s); end;
end
end;
{print info on bottom line}
str (run, s); s := 'run=' + s;
moveto (minX, maxY + 10); outtext (s);
str (batch, s); s := 'batch=' + s + '; SPACE to continue, RETURN to stop';
moveto ((maxX + minX - 8*length (s)) div 2, maxY+10); outtext (s);
str (run+maxX-minX, s); s := 'run=' + s;
moveto (maxX - 8 * length (s), maxY + 10); outtext (s)
end;
procedure plotdot (x, xmin, xmax, w: real; c: integer);
{plot variable x scaled between xmin and xmax in window w, color c}
begin
if x > xmax then x := xmax else if x < xmin then x := xmin;
putpixel (t, round(0.2 * maxY * (w - (x - xmin)/(xmax - xmin))), c)
end;
procedure plot_vars; {plot all interesting variables}
begin
{window 1}
plotdot (x + sqrt (pxx), xmin, xmax, 1, lightgray);
plotdot (x - sqrt (pxx), xmin, xmax, 1, lightgray);
plotdot (x, xmin, xmax, 1, white); {model's x and 1-sigma bounds}
plotdot (xt, xmin, xmax, 1, yellow); {system's x}
plotdot (xopt, xmin, xmax, 1, green); {setpoint}
plotdot (y, xmin, xmax, 1, lightblue);{observation}
{window 2}
plotdot (a + sqrt (paa), amin, amax, 2, lightgray);
plotdot (a - sqrt (paa), amin, amax, 2, lightgray);
plotdot (a, amin, amax, 2, white); {model's a and 1-sigma bounds}
plotdot (at, amin, amax, 2, yellow); {system's a}
{window 3}
plotdot (b + sqrt (pbb), bmin, bmax, 3, lightgray);
plotdot (b - sqrt (pbb), bmin, bmax, 3, lightgray);
plotdot (b, bmin, bmax, 3, white); {model's b and 1-sigma bounds}
plotdot (bt, bmin, bmax, 3, yellow); {system's b}
{window 4}
plotdot (c + sqrt (pcc), cmin, cmax, 4, lightgray);
plotdot (c - sqrt (pcc), cmin, cmax, 4, lightgray);
plotdot (c, cmin, cmax, 4, white); {model's c and 1-sigma bounds}
plotdot (ct, cmin, cmax, 4, yellow); {system's c}
{window 5}
plotdot (u, umin, umax, 5, white) {control}
end;
procedure exit_graphics; {close graphics}
begin
CloseGraph
end;
{extended Kalman filter procedures}
{procedures that describe the behavior of the system}
procedure init_system; {initialize system's parameters}
{these variables describe the system's behavior}
begin
xt := 12.0; {xt := ct + at * xt + bt * u + normal (ft)}
ct := 1.0;
at := 0.9;
bt := 0.1;
ft := fto; {system noise amplitude; pff should agree}
pff := pffo;
vt := vto; {observation noise amplitude; pvv should agree}
pvv := pvvo;
u := 0.0 {action}
end;
procedure step_system; {execute one step of system}
{instead of normal (ft) try also ft * sin (run/5)}
begin
xt := ct + at * xt + bt * u {no unmodelled dynamics}
end;
{procedures that describe the observation}
procedure observe; {do one observation}
begin
if run < 6000 then
y := xt + vt * (sin (run / 5) + sin (run / 155))
else
begin
y := 10.0;
pvv := 1000000
end
end;
{procedures that update the model}
procedure init_model; {initialize model variables}
begin
{parameter estimates}
x := 11.5; {initialize these incorrectly}
c := 1.5;
b := 0.15;
pxx := 4.0; {... but not too incorrectly}
pcc := 4.0; {if you do, see what happens!}
paa := 1.0; {correct: initial paa > (at-a)^2, etc.}
pbb := 1.0;
pax := 0.0; pxa := pax;
pbx := 0.0; pxb := pbx;
pcx := 0.0; pxc := pcx;
pab := 0.0; pba := pab;
pac := 0.0; pca := pac;
pbc := 0.0; pcb := pbc;
qaa := 0.1; {if a changes, by how much?}
qbb := 0.1;
qcc := 0.1;
end;
procedure step_model;
{compute the predicted value of the parameter vector and
its covariance matrix at the next sample interval}
begin
pxx := pcx + a * pxx + x * pax + u * pbx;
pxa := pca + a * pax + x * paa + u * pba;
pxb := pcb + a * pbx + x * pab + u * pbb;
pxc := pcc + a * pcx + x * pac + u * pbc;
pxx := pxc + a * pxx + x * pxa + u * pxb
+ paa * pxx + pax * pax + pff;
x := c + a * x + pax + b * u;
pax := pxa;
pbx := pxb;
pcx := pxc
end;
procedure process_observation; {process one observation y}
{combine prediction and observation into one probability function}
var
t: real;
begin
t := (y - x) / (pxx + pvv);
x := x + pxx * t;
c := c + pcx * t;
b := b + pbx * t;
t := 1.0 / (pxx + pvv);
pxa := pxa - pxa * pxx * t;
pxb := pxb - pxb * pxx * t;
pxc := pxc - pxc * pxx * t;
pxx := pxx - pxx * pxx * t;
paa := paa - pax * pax * t;
pab := pab - pax * pbx * t;
pac := pac - pax * pcx * t;
pbb := pbb - pbx * pbx * t;
pbc := pbc - pbx * pcx * t;
pcc := pcc - pcx * pcx * t;
pax := pxa;
pbx := pxb;
pba := pab;
pcx := pxc;
pca := pac;
pcb := pbc;
end;
{procedures that define how to control}
procedure setxopt; {define setpoint}
begin
xopt := 10.0 + sin (run / 25.0)
end;
procedure control; {compute control; could be improved considerably;
{note that the computation of the control u is based entirely upon
the model's parameters only; NOT on those of the system; NOT on the
observation y}
begin
if b <> 0.0 then
u := (xopt - c - a * x - pax) / b
else
u := 0.0 {if b = 0 then no control is possible, so don't}
end;
{main program}
begin
init_system; {initialize system's parameters}
init_model; {initialize model's parameters}
var_scales; {set bounds for plotting in windows}
run := 0; batch := 0; {startup}
repeat
inc (batch); {next batch of runs}
t := minX - 1; {t = screen x-position for plotting}
init_graphics; {setup axes etc. on screen}
setcolor (red); {outliers signalled by red vertical line}
repeat
inc (run); inc (t); {next run}
step_model; {predict across one sample}
step_system; {... the real thing also}
observe; {observe the system}
process_observation; {collapse of the wave function}
setxopt; {(re)define setpoint}
control; {compute control}
plot_vars; {plot interesting variables}
until t = maxX; {display screen full ?}
repeat until keypressed;
until readkey = chr (13); {stop when RETURN pressed}
exit_graphics {shut down graphics system}
end.
···
a := 0.8;
a := a + pax * t;