The human eye position control system in a rehabilitation setting

Advances in Electrical and Electronic Engineering 118 THE HUMAN EYE POSITION CONTROL SYSTEM IN A REHABILITATION SETTING Y. Nolan a,b) , E. Burke a,b) , C. Boylan a,b) , A. de Paor a,b) a) National University of Ireland, Dublin, Department of Electronic and Electrical Engineering, Belfield, Dublin 4, Ireland b) National Rehabilitation Hospital, Dún Laoghaire, County Dublin, Ireland E-mail: [email protected], Phone: +353 1 7161910, Fax: +353 1 2830921 Summary Our work at Ireland’s National Rehabilitation Hospital involves designing communication systems for people suffering from profound physical disabilities. One such system uses the electro-oculogram, which is an (x,y) system of voltages picked up by pairs of electrodes placed, respectively, above and below and on either side of the eyes. The eyeball has a dc polarisation between cornea and back, arising from the photoreceptor rods and cones in the retina. As the eye rotates, the varying voltages projected onto the electrodes drive a cursor over a mimic keyboard on a computer screen. Symbols are selected with a switching action derived, for example, from a blink. Experience in using this mode of communication has given us limited facilities to study the eye position control system. We present here a resulting new feedback model for rotation in either the vertical or the horizontal plane, which involves the eyeball controlled by an agonist-antagonist muscle pair, modelled by a single equivalent bidirectional muscle with torque falling off linearly with angular velocity. We have incorporated muscle spindles and have tuned them by pole assignment associated with an optimum stability criterion. The dynamics also indicate an integral controller taking its input from a bang-bang element with dead zone. There is, in addition, a pure time delay element involved. Describing Function analysis and simulation demonstrate that in this application the time delay is outside the feedback loop, and is probably associated with set-point generation at a higher level in the brain’s hierarchy of control systems. A second input could be involved at the spindle level, active when tracking predictable target motions. 1. INTRODUCTION It has been known since the mid 19th century that the eye has a standing potential across it (see, for example, the summary by Geddes and Baker [1]) whose sign is dependent on the species. In humans, the cornea is positive with respect to the rear of the eyeball. As the eye moves in its socket, the voltage picked up by a pair of electrodes placed horizontally—on either side of the eye—or vertically—above and below the eye—can be modelled by the simple equivalent circuit shown on Fig. 1, with a typical voltage vs. angle relation graphed on Fig. 2. These have been adapted from Geddes and Baker [1]. We have harnessed these horizontal and vertical electro-oculograms (EOGs) in a communication system for very profoundly disabled people at Ireland’s National Rehabilitation Hospital. The principle of this is illustrated on Fig. 3. Each square on an alphabet board on a computer screen can be highlighted under control of the EOG voltages. During a blink, the vertical EOG voltage is markedly reduced. This is detected and used to write the currently highlighted letter into a message space. Multiple lines of text are accommodated by the “Enter” command. Having harnessed the EOG voltages, we were able to do a very simple experiment and from it deduce a tentative new model for the eye gaze control system in this context. We set up a square on the computer screen, 1.8cm on the side. This was divided into four sub-squares, top left and bottom right black, top right and bottom left white. The subject was asked to focus on the centre of the square, which was then suddenly translated horizontally by an optical angle of 15 º (equivalent to 0.2618 radians) and the resulting EOG recorded. A recording typical of many is shown on Fig. 4. The most significant features of this are the pure time delay of about 0.235 seconds, and the almost perfectly linear excursion over most of the range. This response is almost identical with one given as Fig. 164.5 on page 2481 of Bronzino [2]. Drawing on and extending the work of Stark [3], we propose to explore the control system structure shown on Fig. 5. For each plane of motion, Fig. 1. Genesis of the EOG. The human eye position control system in a rehabilitation setting 119 horizontal and vertical, there is an agonist-antagonist pair of extraocular muscles. Contraction of one rotates the eyeball in the positive θ direction, and of the other in the negative θ direction. We condense these two muscles into a single equivalent muscle, which can rotate the eyeball in either direction. There is no doubt from the physiology of the eye [Davson, 4] that the eyeball muscle torque is controlled by a muscle spindle, although we have not seen that feature invoked in the work of Stark [3] or later authors. However, we quote Davson [4]: “Thus, the spindles in the extraocular muscles are exactly similar to those found in limb muscles and it is therefore impossible to ignore their role in adjusting the force of contraction through a feedback mechanism that indicates the length of the muscle at any moment.” The spindle feedback mechanism is represented by the inner loop in our model. One of our primary objectives is to explore the application of a principle of optimum stability, which we tentatively invoke to tune this inner loop. 2. TUNING THE INNER FEEDBACK LOOP The muscle spindle essentially senses the error, l e, between a locally generated reference value for θ and the current value of θ , and uses it to generate the gross rotational torque, g T, on the eyeball. Drawing on Stark’s [3] work on control of the hand, we take the form of the transfer function relating )(sE l to )(sT g to be 01 2 01 1 )( hshs fsf sC ++ + = (1) We have to “tune” this controller by proposing values for its parameters. As the eyeball rotates (or as a muscle contracts in general), the gross torque falls off with the velocity of contraction. This is generally in a nonlinear fashion, usually taken to be quadratic. However, in the interests of getting a linear dynamical model—at least at this initial stage of analysis—we assume that the fall-off in torque is linear in dtd θ , thus giving the net torque applied to the eyeball as dt d fTT gn θ −= (2) Westheimer [5] has proposed the following model relating θ to n T: θ θθ K dt d FT dt d J n −−= 2 2 (3) Combining this with eqn. (2) and taking the Laplace transform with zero initial conditions gives J K s J Ff s J sT s sG g + + + == 2 1 )( )( )( θ (4) Fig. 2. Approximate dependence of the EOG on θ . Fig. 3. Illustration of a communication system based on the EOG. Fig. 4. EOG in response to a target suddenly translated in the horizontal plane. Advances in Electrical and Electronic Engineering 120 Westheimer [5] gives the values, in SI units, 2 3 12014400 168 102.2 == = ×= − J K J F J (5) With 0=f, these values make )(sG a second order linear system with damping ratio 7.0 = ζ and undamped natural frequency 120 radians per second. We can obtain no guidance from the literature as to what value to assign to Jf , so we propose tentatively but conveniently that it increases the damping ratio to 1 = ζ (critical damping). This gives 72=Jf , and results in 22 )120( 55.454 )120( 1 )( + = + =ss J sG (6) The transfer function relating θ to the local reference value, rl θ , then becomes ( ) ( ) 01 2 01 2 01 1 )120).(( 1 )( )( )( )( fsf J shshs fsf J sG sG s s i i rl +++++ + = = θ θ (7) The denominator here is the characteristic polynomial of the inner loop, whose roots determine the nature of its dynamics. We assume that the function of the spindle is to speed up the eyeball response, so that the roots of the characteristic polynomial lie to the left of the value 120 − = s. We propose to place them all at the location bs 120−= , with 1>b. This gives )( 1 )120).(()120( 01 2 01 24 fsf J shshsbs ++ +++=+ (8) As shown by Cogan and de Paor [6], assigning all the roots of a characteristic polynomial to the same location—which is an extension of the idea of critical damping for a second order system— has an interesting optimum stability property. If all controller parameters but one are held at their nominal values then, as that one is varied through its nominal value, the right-most root is as deep in the left half plane as possible. Dividing Eqn. (8) through by ( ) 2 120+s leads simply to the expressions 424 0 323 1 22 0 1 )120)(386( )120)(412124( )120)(386( )120)(24( −+−= −+−= +−= −= bbbJf bbbJf bbh bh (9) We have explored various values of b and found that 2=b gives a very good final match to the response of Fig. 4. The resulting expression for )( 1 sC is 158400720 22809604.15206 )( 2 1 ++ + =ss s sC (10) For a physical reason, however, we have not used )( 1 sC as given by Eqn. (10), but have developed a close approximation to it. The reason is that Stark’s [3] work strongly suggests that both poles of )( 1 sC should be real, whereas those of )( 1 sC as computed are complex. To retain reality of the poles, we first of all approximated the denominator of )( 1 sC by 22 )400(160000800 +=++ sss . We then scaled 0 f so that the static gain of )( 1 sC , i.e., )0( 1 C, was preserved, and finally scaled 1 f so that 7204.15206800 1 =f. This gave the expression for )( 1 sC actually employed: 160000800 230400016896 )( 2 1 ++ + =ss s sC (11) To give an idea of the accuracy of this approximation, the unit step responses of the two controllers described by Eqns. (10) and (11) are θ DELAY NL C 2 (s) C 1 (s) G(s) nonlinearity controller possible signal from the brain muscle spindle eyeball dynamics θref, from brain + _ + _ e y θ rl e l Fig. 5. Our proposed feedback control system. T g The human eye position control system in a rehabilitation setting 121 shown on Fig. 6, and their Bode magnitude diagrams are compared on Fig. 7. 3. TUNING THE OUTER LOOP Taking guidance from Stark [3], but not following the detailed structure of his model, we now propose that in the forward path of the outer loop we have a nonlinear function followed by a controller. Stark’s [3] treatment suggests that the controller may be modelled by a pure integrator, with transfer function s k sC i =)( 2 (12) However, we differ from Stark in the form of the nonlinear function )(eNL . Stark indicates that this has a dead zone, i.e., 0)( =eNL for 03.0|| <e radians. We preserve this feature, since it represents the fact that no correction need be applied if the image remains focussed onto the fovea, which is the region of greatest sensitivity centred on the optical axis. However, for errors outside this region, Stark has a linear variation of )(eNL with e , and with this feature we have not found it possible to reproduce the linear transition of the response on Fig. 4. With )( 2 sC an integrator, it is hard to escape the conclusion that )(eNL is saturated for 03.0|| >e, and that the input to the inner loop is ramping up linearly during the transition. We consequently adopt the behaviour |||)(| SeNL = for 03.0|| >e, with S a constant, positive value (the saturation level). It is clear that during the transition depicted on Fig. 4, the rate of change of θ is governed by the product i Sk , so that one of these parameters can be scaled arbitrarily. We normalise S to the value 1 = S. With S normalised to unity, only the parameter i k remains to be tuned to get the response closest to that shown on Fig. 4. By simulating the complete system using the Swedish package SIMNON, we have homed in on the value 19= i k (13) With regard to the pure time delay of 0.235 seconds in the response, this has been placed outside the feedback loop, for a reason to be explained in the following section. The resulting response, mimicking conditions in the experiment graphed on Fig. 4, is shown on Fig. 8. 4. OSCILLATIONS IN THE SYSTEM In seeking to place the time delay correctly, we noticed that even a very small delay (to be quantified below) placed within the loop led to continuous oscillations. The same happens if i k is raised sufficiently, or if the width of the dead zone in NL is decreased sufficiently. In all these cases the output oscillation is almost sinusoidal, and this suggests that the classic Describing Function technique [Elgerd, 7] can be used to explore the conditions for oscillation. In applying the describing function, we assume that the error signal, e , in the outer loop is oscillating (almost) sinusoidally: )sin( tMe ω = (14) For D M ≤ , where D is the half-width of the dead zone ( 03.0 = D in this case), 0 = y, which is of no interest. However, for D M > , y is periodic, of period ωπ 2=T, but non-sinusoidal. In fact it is a Fig. 6. Unit step response of the two spindle transfer functions: Eqn. 10 (a) and Eqn. 11 (b). Fig. 7. Bode magnitude diagram of the two spindle transfer functions: Eqn. 10 (a), Eqn. 11 (b). Fig. 8. Overall step response of the complete system. Advances in Electrical and Electronic Engineering 122 Fig. 9. Frequency response locus of Gf(s). gapped square wave. If we expand y in a Fourier series, we find that, due to the symmetry of NL , it has only odd harmonics, and all of those above the fundamental frequency are effectively filtered out in passing through the integrator and the inner feedback loop. Thus, as far as its effect on the observed oscillation is concerned, we can approximate y by its fundamental component, readily calculated as )sin(. 4 22 t M DMS y f ω π − = (15) Comparing Eqns. (14) and (15), we see that for this almost sinusoidal oscillation, NL can be characterised by a gain, referred to as its Describing Function, DF : DM M DMS DF > − =,. 4 2 22 π (16) Plotted as a function of M , this reaches a peak value of D S DF π 2 max = (17) for 2DM = (18) In the extreme case 221.21 max == DFDF (for 1 = S, 03.0 = D), we can use an informal application of Nyquist analysis [Elgerd, 7] to explore conditions for oscillation. To do this we plot the frequency response locus of the function )( 19 )( max sG s DFsG if       = (19) This is shown on Fig. 9, which also indicates certain “robustness margins”. These tell by how much the frequency response locus should be modified in various ways to pass it through the critical point )0,1(−, which corresponds to sustained oscillations. Thus, for example, the gain margin of 1.888 means that if i k is multiplied by this factor, oscillation sets in. This is illustrated on Fig. 10. Similarly, the delay margin of 0.0047 means that if a pure time delay of this amount is placed in the forward path of the outer loop, oscillation ensues at the angular frequency 108.3 radians per second. The phase margin of 29.11°, also effective at the angular frequency 3.108 = ω radians per second, means that if an extra dynamic element were inserted in the forward path, with a phase lag of 29.11° at the stated angular frequency, the frequency response locus would pass through the critical point. We have confirmed the findings with respect to gain margin and delay margin by simulation. The very small delay margin is the reason why we had to place the observed time delay of 0.235 seconds outside the loop. 5. DISCUSSION As an offshoot to development of a communication system for a disabled people, based on the EOG, we performed a simple experiment to track a suddenly translated target, and used it to produce a new feedback model for the eye gaze control system in this application. Muscle spindles have been incorporated into this type of model for the first time. In the absence of any detailed information, the resulting inner feedback loop has been tuned by a pole-assigning procedure, associated with a principle of optimum stability. In the outer loop, the nonlinearity has been modified from that suggested by Stark [3] and, on tuning the integral controller gain through simulation experiments, a very close match to the experimental result. Robustness margins have been explored with the help of Describing Function analysis, and this has been decisive in placing the time delay outside the control loop. We hope that this exercise in applying several techniques from Control Theory to a living system will be of interest to the Biomedical Engineering community. Fig. 10. Oscillations due to setting ki just above its threshold value. The human eye position control system in a rehabilitation setting 123 Acknowledgement The authors are very grateful to the Irish Research Council for Science, Engineering and Technology and to the National Medical Rehabilitation Trust for their support of this work. REFERENCES [1] L.A. Geddes and L.E. Baker, Principles of Applied Biomedical Instrumentation, John Wiley & Sons, New York, 1989 (3rd Edition). [2] J.D. Bronzino, editor-in-chief, The Biomedical Engineering Handbook, CRC Press, Boca Raton, Florida, 1995. [3] L. Stark, Neurological Control Systems: Studies in Bioengineering, Plenum Press, New York, 1968. [4] H. Davson, The Physiology of the Eye, Little, Brown, Boston, 1963. [5] G. Westheimer, Mechanism of saccadic eye movements, AMA Archives Ophthalmology, vol. 52, 710-724, 1954. [6] B. Cogan and A. de Paor, Optimal stability and minimum complexity as desiderata in feedback control system design, IFAC Conference, Control Systems Design, Bratislava, Slovakia, June 18-20, 2000, pp. 51-53. [7] O.I. Elgerd, Control Systems Theory, McGrawHill, New York, 1967.