UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Modelling neural transdution Dean, Douglas Philip 1984

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-UBC_1985_A1 D42.pdf [ 5.31MB ]
Metadata
JSON: 831-1.0096542.json
JSON-LD: 831-1.0096542-ld.json
RDF/XML (Pretty): 831-1.0096542-rdf.xml
RDF/JSON: 831-1.0096542-rdf.json
Turtle: 831-1.0096542-turtle.txt
N-Triples: 831-1.0096542-rdf-ntriples.txt
Original Record: 831-1.0096542-source.json
Full Text
831-1.0096542-fulltext.txt
Citation
831-1.0096542.ris

Full Text

MODELLING NEURAL TRANSDUCTION by DOUGLAS PHILIP DEAN B.A.SC, Engineering Physics, The University of British Columbia, 1978 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES (Department of Electrical Engineering) We accept this thesis as conforming fce~t.he required standard THE UNIVERSITY OF BRITISH COLUMBIA August 1984 (c) Douglas Philip Dean, 1984 In p r e s e n t i n g t h i s t h e s i s i n p a r t i a l f u l f i l m e n t of the requirements f o r an advanced degree at the U n i v e r s i t y o f B r i t i s h Columbia, I agree t h a t the L i b r a r y s h a l l make i t f r e e l y a v a i l a b l e f o r r e f e r e n c e and study. I f u r t h e r agree t h a t p e r m i s s i o n f o r e x t e n s i v e copying of t h i s t h e s i s f o r s c h o l a r l y purposes may be granted by the head of my department or by h i s or her r e p r e s e n t a t i v e s . I t i s understood t h a t copying or p u b l i c a t i o n of t h i s t h e s i s f o r f i n a n c i a l gain s h a l l not be allowed without my w r i t t e n p e r m i s s i o n . Department of •g.uecm^OcC & J 6 3 > ) U ^ a £ ^ »-JCq> The U n i v e r s i t y of B r i t i s h Columbia 1956 Main Mall Vancouver, Canada V6T 1Y3 Date 1 3 L - / 2 J ABSTRACT Neural transduction i s the process by which neurons convert externally applied electrical energy into stable, propagating voltage pulses. Given some stimulus waveform with particular parameters such as duration, phase delay, etc., there i s a minimum stimulus amplitude required in order for transduction of the waveform to result in an active neural response. The minimum amplitude for excitation, the threshold amplitude, i s a strong func-tion of many variables including stimulus waveshape, frequency and duration. This study reveals some details of the threshold characteristics of the Frankenhaeuser-Huxley (FH) model of myelinated nerve. These threshold trans-duction characteristics are studied with the aid of phase-space analysis, and are used to produce a model of neural excitation which i s c l i n i c a l l y applic-able to human nerve. The f u l l FH system of equations i s used to predict threshold behaviour for in-vivo human nerve, and the predictions are shown to be in good agreement with c l i n i c a l l y obtained threshold data. The study concludes by examining some of the additions to the FH model which would be necessary to model the accumulation of extra-nodal potassium ions. i i TABLE OF CONTENTS Page ABSTRACT i i TABLE OF CONTENTS i i i LIST OF FIGURES vi ACKNOWLEDGEMENT v i i i PREFACE ix 1. INTRODUCTION 1 1.1 Modelling Neural Transduction 2 1.2 Definition of Terms 4 1.3 Problem Definition 5 1.4 The Basis for FES as a Successful Mode of Treatment 7 1.5 Nerve -vs- Muscle Stimulation in FES 8 1.6 Factors Influencing Stimulation 9 1.6.1 External Factors: Stimulus Parameters, Stimulus Design 9 1.6.2 Electrode Tissue Interface 11 1.6.3 Internal Factors: Membrane Transduction 12 1.6.4 Internal Factors: Motor Unit Responses 12 1.7 Choosing a Computational Model of Nerve 13 2. PHASE-SPACE ANALYSIS OF THE FRANKENHAEUSER-HUXLEY MODEL 17 2.1 Introduction 17 2.2 Numerical Solution of the FH System 19 2.3 Classification of the FH System 21 2.4 Phase Space Analysis of the FH System 22 2.4.1 Choice of Phase Space - A Reduced System 23 2.4.2 Isoclines and Cr i t i c a l Points in the Vm Phase Plane . 25 2.4.3 Threshold Characteristics of the FH System 28 2.4.4 Phase Plane Variations with n, p, h 30 2.4.4.1 Potassium Activation Variable n 31 2.4.4.2 Unspecific Current Activation Variable p . . 31 i i i Page 2.4.4.3 Sodium Inactivation Variable h 31 2.4.4.4 Sodium and Potassium Permeability 35 2.5 Determining Threshold Stimulus Parameters 36 2.5.1 The Separatrix Definition of Threshold 42 2.6 Discussion 47 2.7 Summary 50 3. OPTIMIZATION OF NEURAL STIMULI BASED UPON A VARIABLE THRESHOLD POTENTIAL 51 3.1 Introduction 51 3.2 The Classical Strength-Duration Curve 53 3.2.1 Limitations of the RC Strength-Duration Curve . . . . 55 3.3 The Threshold Potential Function 57 3.3.1 Monophasic Stimuli 58 3.3.2 Biphasic Stimuli 58 3.4 Threshold Current and Charge 60 3.5 Stimulation with Minimum Power 61 3.6 Stimulation with Minimum Charge 66 3.7 Clinical Determination of Threshold Parameters 68 3.8 Discussion 71 3.9 Conclusions 74 4. SIMULATION OF HUMAN MEDIAN NERVE THRESHOLD RESPONSE 75 4.1 Introduction 75 4.2 Enhancement and Extension of the Model 77 4.2.1 Temperature 78 4.2.2 Passive Tissues and Adjacent Nodes 81 4.3 Experimental Design 86 4.4. Experiments 89 4.4.1 Percutaneous Stimulation 91 4.4.1.1 Selection of Model 91 4.4.1.2 Electrophysiological Techniques 91 4.4.1.3 Stimulus Waveforms 93 iv Page 4.4.1.4 Reference Selection 93 4.4.1.5 Procedure 94 4.4.1.6 Comparison with the Model 95 4.4.2 Transcutaneous Stimulation 98 4.4.2.1 Selection of Model 98 4.4.2.2 Methods and Procedures 100 4.4.2.3 Stimulation 100 4.4.2.4 Comparison with the Model 101 4.5 Discussion 103 4.6 Conclusions 107 5. SIMULATION OF NODAL POTASSIUM ACCUMULATION 109 5.1 Introduction 109 5.2 Models of K+ Accumulation 110 5.2.1 The Three Compartment Model 110 5.2.2 Diffusion in an Unstirred Layer Model I l l 5.3 Choice of a IT*" Accumulation Model for the FH System I l l 5.4 The Affects of K+ DUSL Accumulation on the FH System . . . . 112 5.5 The Influence of K+ Accumulation on Multiple Pulse Thresholds 117 5.6 Discussion 121 6. SUMMARY 123 REFERENCES 128 APPENDIX I 134 v LIST OF FIGURES Page 1.1 Elements of Neural Stimulation 16 2.1 Standard Response of the FH Model 20 2.2 Quasi-Threshold Configuration of the FH Phase Plane 24 2.3 Resting State Isoclines of the Vm Phase Plane 27 2.4 Resting State Phase Plane Detail 29 2.5 Affects of Changing 'n' on the Vm Phase Plane 32 2.6 Affects of Changing 'p' on the Vm Phase Plane 33 2.7 Affects of Changing 'h' on the Vm Phase Plane 34 2.8 Vm Phase Plane at Peak P„ and Peak P., 37 Na K 2.9 Stimulus-Response Curves for 120 us Monophasic Square Stimuli at 20°C and 15°C 39 2.10 Subthreshold Response to 5 us Biphasic Square Stimulus 0.0 us Phase Delay 41 2.11 Threshold Amplitudes of Biphasic Stimuli Determined by the Separatrix Definition 44 2.12 Comparison of Threshold Calculation Times for Separatrix and Level Definitions 46 2.13 Stimulus-Response Curves for 120 us Monophasic Square Stimuli at 2.5°C, 20°C, and 40°C 49 3.1 RC Model Strength-Duration Curves Compared to FH Model Threshold Amplitudes 56 3.2 Threshold Potentials for Biphasic Square Stimuli from the FH Model 59 3.3 RC-Variable Threshold Potential Strength-Duration Curves Compared to FH Model Threshold Amplitudes 62 3.4 Pulse Widths Minimizing the Normalized Damage Function 69 v i Page 4.1 Equivalent Circuit of the Three-Node FH Model 82 4.2 Distortion of Injected Current Waveforms . . . 84 4.3 DISA 15E07 Stimulus Waveform 85 4.4 Near Nerve Recordings 90 4.5 Assumed Geometry for Percutaneously Applied Stimuli 92 4.6 Model and Clinical Thresholds for Percutaneous Stimuli 96 4.7 Assumed Geometry for Transcutaneoujsly Applied Stimuli 99 4.8 Model and Clinical Thresholds for Transcutaneously Applied Stimuli 102 4.9 Cl i n i c a l Thresholds for Percutaneously and Transcutaneously Applied Stimuli 106 5.1 Potassium Current Flow in Response to Standard Stimulus . . . . 114 5.2 FH Variables Affects by Potassium Ion Accumulation . 115 5.3 Phase Plane Changes due to Potassium Ion Accumulation 116 5.4 Phase Plane Detail 118 5.5 Threshold Delay Times 120 v i i ACKNOWLEDGEMENTS I wish to thank Peter Lawrence for his assistance and guidance during the course of this study. Without his encouragement this thesis would not have been finished. I must also thank Theresa Barber and Dave MacDonald for the part they played in the completion of the work described on these pages. In addition, Paul Kendrick and Steve Mennie spurred my creative urges with observations and comments from a non-technical point of view — their input was more important than they realize. A special acknowledgement is due to Miss Heart, who really started a l l of this in 1962. Finally, many thanks to B i l l March for advice, stern warnings, and much much more. v i i i PREFACE This thesis i s presented as four self-contained chapters. The presenta-tion of each chapter is progressively less theoretical and more c l i n i c a l , beginning with a mathematical analysis of the properties of the Frankenhaeuser-Huxley system of equations, and leading to the presentation of c l i n i c a l data demonstrating the applicability of model results to human nerve stimulation. With the exception of the fourth chapter and the Introduction, a l l chap-ters have either appeared as papers in the IEEE Transactions on Biomedical Engineering, or are awaiting publication, or have been submitted for publica-tion in that journal. Consequently, each chapter contains i t s own introduc-tion, summary, conclusions and discussion, and may be read without reference to preceeding or following chapters. ix 1 1. INTRODUCTION The motivation for this work has been a perceived need to understand more about the threshold characteristics of myelinated human peripheral nerve fibres. Understanding the nature of neural threshold i s an important goal since many significant applications of neural control (as outlined below) implicitly involve making decisions regarding stimulus waveshapes and parameter settings which are ultimately predicted upon aspects of neural threshold. It i s d i f f i c u l t to perform experiments on the human peripheral nervous system. The system i s delicate, important, and to a large extent nonregener-able. Hence the primary source of investigation in this work was not in-vivo or in-vitro human nerve, but rather a computational model implemented on a di g i t a l computer. By using such an 'electronic' nerve, i t was possible to perform some threshold 'experiments', or experimental analogues, which would not have been possible to perform on a human subject. Experiments with the human nervous system were undertaken only to obtain data which could be used to provide an indication of the relevance of the computational model to the system being simulated. No experiments were performed on human nerve for the purpose of gathering physiological data which would be used to alter model parameters. The question of model applicability to human nerve i s , of course, c r i t i -cal since no matter how convenient the computational model may be, i t s results are irrelevant i f one wishes them to provide insights into the behaviour of human nerve, but one finds that model dynamics have no resem-2 blance to those of human peripheral nerve. The question of choosing a neural model which has the greatest probable applicability to the human case is addressed at the end of this introduction. A transducer i s a device which converts the form of some applied signal into an alternate, usually e l e c t r i c , signal. The process by which the input signal i s converted into the output signal i s known as transduction. The neurons which compose the peripheral nervous system in man behave like trans-ducers when subjected to externally applied electrical energy; however, the neuron i s a rather unusual transducer of external electrical signals in that i t operates on an ' a l l or none' protocol. If the input signal (the external-ly applied el e c t r i c a l energy) meets certain c r i t e r i a , the neuron-transducer produces an output signal which has a characteristic amplitude and shape, and which w i l l travel from the site of transduction at a characteristic velocity. This stable, propagating voltage pulse i s known as the action potential. If the input signal does not meet the threshold c r i t e r i a , no stable propagating action potential w i l l be produced. 1.1 Modelling Neural Transduction Regardless of the type of external neural stimulation being applied, the properties of the neural membrane which f a c i l i t a t e action potential genera-tion are of c r i t i c a l importance. The purpose of this work i s to analyse the threshold transductive properties of the myelinated axon with the aid of the Frankenhaeuser-Huxley model [1], and to discuss some of the resulting c l i n i -cal implications. 3 Within the peripheral nervous system, neurons (nerves) are arranged in bundles or f a s i c u l i which are surrounded by connective tissue sheaths. Thousands of neurons, or neural elements, may be present i n a single neural bundle. The neuron is composed of a nerve c e l l body and i t s processes which are called dendrites and axons. The axonal process serves to conduct impulses to other neurons or to muscle fibres and gland c e l l s , while the dendritic process serves as a point of termination for axons of other neurons. The reader i s referred to an introductory textbook on the human nervous system such as that by Barr and Kiernan [2] for further details. A block diagram, Fig. 1.1, shows how the properties of the neural ele-ment are related to other factors influencing external neural stimulation. In this diagram, I refers to the externally applied current, i n units of t ll current density, which serves as an input signal to the n neuron; and t h i s the transduced response of the n neuron to 1^. Although Sr is s t r i c t l y a voltage signal, i t i s considered here as a dimensionless binary signal with the two binary states representing the a l l or none condition of transmission. A stimulus current IQ at some potential Vq i s applied through a set of stimu-lus electrodes. Assuming the use of a constant current stimulator, the impedance of the electrode-tissue interface (ETI) does not affect the stimu-lus current, and so a current of I passes through the electrodes into surrounding tissues. The action of this current on the ETI can produce stimulation by-products which can, in turn, influence both ETI impedance and the properties of the tissues near the electrodes. As the stimulus current passes through passive tissues surrounding the neural elements, some frac-4 t i o n , T-n» i s diverted to enter each of n myelinated neural elements through a node of Ranvier. If the current I entering the element meets neces-sary threshold requirements, then the neural element w i l l be induced to generate an active response - the action potential. This active response, S^, w i l l propagate from the site of activation with a characteristic velocity and shape and w i l l be processed in the central nervous system i f the neural element was responsible for sensory input, or Sr may result in activation of a motor unit i f the neural element was a moto-neuron. 1.2 Definition of Terms Electrical stimulation may be applied to the body either for therapeutic or diagnostic purposes as pointed out by Roberts e t . a l . [3]. The application of el e c t r i c a l energy to human neurological tissues for the purpose of evoking active neural responses may be called applied neural control or neuro- modulation. The therapeutic goals of applied neural control, the means of applying and delivering the electrical energy to the body, and the targetted tissues, provide convenient means of classifying the various neuromodulation scenarios. The two most general classes of neuromodulation are created by two dis-tinct classes of therapeutic application - one where the intention i s to enhance weak neurological signals or generate non-existent neurological signals - the other where the purpose i s to mask or supress undesired neural a c t i v i t y . 5 In cases where disease or trauma have damaged c r i t i c a l neural pathways, the resulting sensory or motor deficits may be somewhat alleviated by applied neural control. Such control or information signal enhancement has been called neuroaugmentation and functional el e c t r i c a l stimulation. Undesired neural activity with few exceptions means perception of pain, either chronic or acute. In these cases the object of neuromodulation is to eliminate the perception of pain. The mechanism of pain suppression by neural stimulation i s not well understood, but i t i s thought that the release of endogenous opiates may be involved. Further sub-classifications of therapies in the context of the above two classes of neuromodulation can be made on the basis of the method of applying the el e c t r i c a l stimulus. In spinal cord stimulation, electrodes are surgi-cally implanted to stimulate the dorsal column of the spine, similarly cortical stimulation refers to stimulation of the brain via surgically inserted electrodes. Percutaneous stimulation, on the other hand, involves the use of needle-like electrodes inserted through the skin, without surgical intervention, for the purpose of stimulating peripheral nerves, while trans- cutaneous stimulation stimulates nerves through unbroken skin. 1.3 Problem Definition Considering only functional electrical stimulation (FES), and not stimu-lation for pain r e l i e f , the goal of neuromodulation may be said to be the generation of functional responses in 'paralysed' end organs which are indis-tinguishable from the responses of identical 'unparalysed' organs, through application of externally generated control signals. 6 Human nerves can be induced to transmit control signals by application of externally generated stimuli. Fortuitously, such control signals are identical to those which are naturally generated, and hence the organs under control of the a r t i f i c i a l l y stimulted nerves w i l l respond in an identical manner to the a r t i f i c i a l control signals. In principle then, one should be able to restore lost motor or perceptual function by stimulation of the nerves controlling the organs in question. In practice, such a restoration of lost function embodies a large number of c r i t i c a l factors from stimulator design and electrode positioning, to the transductive properties of the neural membranes being stimulated. In this thesis, we study the c r i t e r i a which must be met by the stimulus current I passing into the nC^ neural element so that the neural membrane gives an active response. How the resulting signals are processed, or how they activate motor end units are not studied. The properties of surrounding passive tissues are addressed only insofar as is required to examine the correspondence of c l i n i c a l data to model predictions. It is not the purpose of this work to model in detail passive tissues, electrode-tissue interfaces or geometries of stimulating electrodes. Rather, the threshold properties of neural elements in response to external stimuli and some of the resulting c l i n i c a l implications are investigated. The FES c r i t e r i a may be stated as follows: 1. Nerves must be stimulated by pulses generated from a stimulator which meets established safety requirements; 7 and which produces comfortable, effective trains of electrical stimuli. 2. The stimuli must be applied through electrodes which minimize discomfort or damage to surrounding tissues and which have a geometry which is optimal for the effect desired. 3. The stimulus parameters selected (frequency, duration, etc.) must simultaneously activate the target neurons, and be perceived as comfortable by the patient. 4. The induced control signals must have the spatial and temporal properties to produce smooth regular contraction of the muscle which is innervated by the neural elements being stimulated; or in the case of sensory organs, that an appropriate sensory response i s generated. In practice, 4 i s extremely d i f f i c u l t to achieve to the extent that there i s no noticeable difference between a naturally generated muscle contraction, and one generated by techniques of functional electrostimulation. 1.4 The Basis for FES as a Successful Mode of Treatment The principle of neural control signal indistinguishability is one of the fundamental axioms of functional electrostimulation. It i s well known that the fundamental unit of neural communication, the action potential, can be generated by subjecting an in-vivo nerve to a changing electric f i e l d . The physiological mechanisms of the neural membrane which are activated by these ' a r t i f i c i a l ' means to produce a stable propagating action potential are 8 Identical to the mechanisms of 'natural' action potential generation and conduction. Away from the immediate v i c i n i t y of the applied electrical stimulus there i s therefore no distinguishable difference between naturally generated and a r t i f i c i a l l y generated action potentials. Hence, any organ under control of the neural element being stimulated w i l l respond to the stimulus induced action potential in exactly the same manner as i t would respond to a similar naturally generated action potential. 1.5 Nerve -vs- Muscle Stimulation in FES When applying FES to patients, i t is possible to produce the desired functional effects by direct stimulation of nerve, or by stimulation of skeletal muscle, as stated by Mortimer [4]. Both approaches have been used by many different research groups i n the past, but direct nerve stimulation has proven to have several advantages compared to muscle stimulation: a) Lower stimulus amplitudes are required - this prolongs stimulator battery l i f e and reduces possible stimulation side effects. b) Electrodes can be located in areas which do not move -electrode displacement (or breaking in the case of implants) i s common when the electrodes must be located on the body of the muscle which i s being stimulated. c) One muscle function can be controlled by a single nerve. The chief problem of direct nerve stimulation - that of fatigue at the neuro-muscular junction - occurring only in the case of implanted stimulus elec-trodes, has been largely eliminated by the technique of 'karusellstimulation' 9 developed in Austria by Holle et a l . [5]. With this technique, a cyclic inhomogeneous elec t r i c a l f i e l d i s generated by implanted electrodes surround-ing the nerve to be stimulated. The rapid variation of the cyclic f i e l d eliminates neuro-muscular junction fatigue by inducing a more equitable dis-tribution of control signal conduction throughout the target nerve. 1.6 Factors Influencing Stimulation Factors of interest in FES may be roughly organized into three cate-gories [4]. These classifications are defined by factors external to the body - mainly those concerned with stimulator design and function (stimulus waveform, stimulation parameters); factors influencing the electrode-tissue  interface such as electrode impedance and stimulation by-products; and fac-tors wholly internal - motor unit responses to stimuli, transductive proper-ties of the neural membrane, and 'central processing' of signals from afferent sensory neurons. These factors and their Inter-relationships are shown in Fig. 1. Further classifications are usually made based upon the nature of the neural elements being stimulated - central nervous system (CNS) stimulation, or peripheral nervous system (PNS) stimulation. Within the PNS class of applications, there are further sub-classifications according to stimulation of either efferent or afferent pathways. 1.6.1 External Factors: Stimulus Parameters, Stimulator Design It i s accepted by most researchers that charge per stimulus i s a crucial factor in stimulation, thus implying that complete control over stimulus current i s necessary. This further implies that a constant current stimula-10 tor i s preferable to a constant potential stimulator since in the la t t e r , current output i s a strong function of tissue and electrode impedance (Butikofer and Lawrence [6]). Most FES stimulators are thus of the constant  current variety. To avoid dangerous increases in current stimulator driving potential which can occur in the case of an increasing electrode impedance, an over-voltage detection circuit should be incorporated into the stimulator design. These two considerations, constant current stimulation and over-voltage detection, are recognized by nearly a l l researchers to be the two most fundamental requirements of any FES stimulator design. Although stimulus waveform and parameter investigation i s s t i l l an active area of research, certain parameter ranges and stimulus waveforms are agreed upon as most effective. In general, pulsatile biphasic stimuli are preferred, with stimulus durations ranging from 20 [is to 300 us for transcu-taneously applied stimuli, and 1 |is t o 300 us for percutaneously applied stimuli (Vossius [7]). Pulsatile rather than sinusoidal stimuli are pre-ferred since with sinusoidal stimulation charge delivered per stimulus phase is dependent upon stimulus frequency. This i s not the case with pulsatile stimulation. Stimuli should be biphasic and charge balanced in order to avoid perturbing electro-chemical equilibria at the electrode-tissue inter-face (McHardy et a l . [8], Pudenz et a l . [9]). Recent studies have shown that different stimulus waveforms may be arranged in a heirarchy of efficacy with respect to minimal stimulation amplitudes and minimal noxious sensations during stimulation (Janko and Trontelj [10]). Although there i s a general concensus regarding the form of the most favourable stimulus waveshape, there 11 are s t i l l many aspects of stimulus parameter selection which are not well understood and which require further research. 1.6.2 Electrode Tissue Interface Electrode-tissue interface (ETI) considerations have important implica-tions for stimulator design. The principal factors of concern relate to the interface impedance, and to the dissemination of electro-chemical reaction  products to tissues in contact with the stimulating electrodes. If the electrode impedance during stimulation becomes too high for reasons of diminishing electrode contact, local regions of high stimulus current density may occur resulting in tissue damage due to ohmic heating and associated painful sensations (Burton and Maurer [11], Mason and MacKay [12]). Charge imbalanced stimuli can perturb the thermodynamic electro-chemical equilibria of the ETI resulting in corrosion of the stimulus electrode. The same phenomena may also occur with charge balanced stimuli which have exces-sively long phase durations. Transcutaneously, this equilibrium imbalance may manifest i t s e l f as local reddening and i r r i t a t i o n of the skin [12], while in the case of implanted electrodes toxic electro-chemical by-products may induce local tissue breakdown (Pudenz et a l . [13]). Additionally, the different current dispersion influences of various ETI configurations can significantly alter the appropriate magnitudes of stimulus current waveforms, even changing the neural response patterns as functions of stimulation parameters. This factor must be considered when selecting an FES stimulator for a particular application. 12 1.6.3 Internal Factors; Membrane Transduction The threshold properties of a neural element are of prime importance in stimulus parameter selection, stimulus waveform selection and electrode posi-tioning. Ideally one would like to induce neural activity with the minimum possible stimulus amplitude in order to minimize both ETI problems and to minimize current drain on the batteries powering the electrical stimulator. Based upon modelling studies and experimental verification ([14]-[16]), threshold stimulus amplitude l o c i as functions of stimulus waveform parameters are known and can be used to select stimulus parameter ranges to optimize stimulus efficacy (Gorman and Mortimer [17]). In addition, some recent work has detailed the form which an optimal stimulus should take [18]. 1.6.4 Internal Factors; Motor Unit Responses Afferent nerve stimulation, achieved by activation of skin receptors and the largest diameter A-p fibres of the afferent nerve bundles, triggers spinal reflex arcs and may be used to generate complex movements. It i s well known that afferent fibre FES has beneficial secondary effects on the a c t i -vated motor units, decreasing muscle atrophy, increasing both muscle tone and blood flow, as well as decreasing inactivity osteoporosis in bone [2], [7] and [19]. Direct stimulation of efferent nerve fibres produces useful muscle con-tractions, and has been demonstrated with the electronic peroneal brace (Vodovnik et a l . [19]), with gait assists in hemiplegics (Stanic et a l . [10]), and with standing assists for paraplegics using transcutaneously applied FES signals [21], and implanted electrodes (Sthor et a l . [22]). In 13 addition i t has been shown that stimulation of efferent nerve fibres improves volitional control (Gracanin [23]) and reduces spasticity (Rebersek et a l . [24], Bowman and Bajd [25]). A combination of afferent and efferent FES has been applied to hemi-plegics, using efferent FES during the gait swing phase with simultaneous afferent stimulation of the non-moving leg to produce a flexion response [21]. 1.7 Choosing A Computational Model of Nerve A large number of computational models of nerve have been presented in the literature. Cole, [26], provides a comprehensive review of the most important neural models. In general, a model of nerve which accounts for the observed threshold characteristics must consist of a capacitance in parallel with a time-varying non-linear resistance. The Hodgkin-Huxley [2 7] model of unmyelinated nerve was the f i r s t model to separate the membrane resistance into distinct components - one for each ionic species, Na+, K"1" and C l ~ . Perhaps the most well known models of nerve are the afore-mentioned Hodgkin-Huxley (HH) model of squid unmyelinated axon; the Frankenhaeuser-Huxley (FH) model of myelinated nerve from the frog Xenopus leavis [1]; Dodge's model of Rana pipiens myelinated nerve [28]; and Noble's model of mammalian Purkinje fibre [29]. Noble's model is unsuited to investigations of peripheral nerve since i t emulates the dynamics of the pacemaker nerve fibres in mammalian heart whose behaviour is strikingly different from peripheral nerves. Dodge's model is not scaled to unit area and i s , i n fact, a model of a single node of the many 14 which were investigated ("node 7"), and hence is unsuitable for our applica-tion. Although the HH model is valid, strictly speaking, only for the giant axon of squid at 2.5°C, and the FH model is valid only for the nerve fibres of Xenopus laevis, a more optimistic view may be that the HH model is a generic model of unmyelinated nerve, and the FH model one of myelinated nerve. Such a hypothesis cannot be proven, but many examples abound in the literature which seem to indicate that the HH and FH models emulate many of the mechanisms of excitation in species other than those which supplied the data bases for the original models. For example, Fitzhugh [30] showed using dimensional analysis that results from one nerve model can be transformed to apply to another fibre of similar structure having different model constants. The experimentally verified linear diameter-velocity relation of myelinated fibres was derived from the FH equations by Goldman and Albus [31], while Buchthal and Rosenfalck [32] determined that the relation between temperature and conduction velocity in human sensory nerve was similar to that in frog nerve and in the vagus and saphenus nerves of cat. Although encouraging, the above examples by no means prove that such models can be used to produce results which predict human nerve dynamics. Since we are interested in the large diameter myelinated fibres of the human peripheral nervous system, the FH equations were chosen as the basis of the computational model in this study. In chapter 2, the dynamics of the FH model are studied using phase-space analysis. The analysis of chapter 2 leads to a precise definition of threshold which is useful in determining threshold amplitude loci of various stimuli as functions of stimulus parameters. Chapter 3 makes use of the threshold definition in chapter 2 to suggest an algorithm for selecting 'optimal' stimulus waveforms and parameter settings. In chapter 4, we w i l l present c l i n i c a l results from human subjects which compare favourably to analogue FH model computations. These results w i l l provide some j u s t i f i c a t i o n for the decision to employ the Frankenhaeuser-Huxley model as a means of emulating the transductive properties of myelinated human peripheral nerve. 16 CENTRAL PROCESSING (SENSORYI CONSTANT CURRENT STMULATOR ELECTRODE-TISSUE MTERFACE I I I JL. TlSSUE STMULATON BVPROOUCTS NEURAL ELEMENT X MOTOR UNIT ACTIVATION Figure 1.1 Elements of Neural Stimulation 17 2. PHASE SPACE ANALYSIS OF THE FRANKENHAEUSER-HUXLEY MODEL 2.1 INTRODUCTION The Frankenhaeuser-Huxley (FH) equations are a system of five non-linear ordinary differential equations which describe the space-clamped response of a myelinated toad (Xenopus laevis) nerve to current stimuli at a node of Ranvier. The equations were empirically derived by Frankenhaeuser, Dodge and Moore over a number of years using the voltage clamp technique in a series of experiments [33] - [42]. The results of the experiments and the derived set of equations are summarized and discussed by Frankenhaeuser and Huxley in [1], and are repeated in Appendix I. The equations give an expression for the f i r s t derivative with respect to time of transmembrane potential as a sum of ionic current components of sodium, potassium, an unspecific component (primarily carried by sodium ions), a stimulus current and.a leakage current. The primary ionic currents (K+, Na+ and unspecific) are expressed as a product of time-varying perme-a b i l i t y of the membrane to the ionic species in question, and a non-linear Goldman factor. The Goldman factor, which describes the observed rectifying properties of a semi-permeable membrane, can be derived from the assumption of a constant f i e l d membrane. It was observed by Dodge and Frankenhaeuser [33], [34], that peak i n i t i a l current densities caused by step depolariza-tions of the membrane departed significantly from a linear form in the region of the largest depolarizations. In [34] they proposed that the ionic cur-rents would best be described by the constant f i e l d Goldman equation [43]. Introduction of the delayed unspecific ionic current component, and description of the ionic currents using Goldman factors are the fundamental 18 differences between the FH system and the well known Hodgkin-Huxley equations for unmyelinated nerve [2 7]. In the Hodgkin-Huxley system, a linear relation between membrane current and transmembrane potential was used in conjunction with a time varying ionic conductivity. Numerical solutions of the FH equations were used by Butikofer and Lawrence [44] to calculate minimum amplitudes of biphasic square current stimuli required to e l i c i t action potentials in the FH system. In the same paper, the authors demonstrated the a b i l i t y of FH model solutions to predict trends of experimental data corresponding to bulk stimulation of non-space-clamped human nerves. In a later paper [6] the same authors used their previously calculated threshold data to determine the thermal energy generated in the skin during the course of a transcutaneously applied threshold stimulus. The FH model was used by McNeal [16] to investigate current stimulation of a multi-node myelinated nerve. In order to reduce the complexity of his multi-node model, McNeal used a simplified passive form of the FH equations and calculated membrane current and potential for small depolarizations (<15 mV) where the assumption of passive RC response by the FH system i s v a l i d . Mortimer [4] employed the model described by McNeal [16] to determine the threshold dependence of myelinated nerve as a function of axon diameter and position with respect to stimulating electrode. The threshold character-i s t i c s were investigated for monophasic constant-current stimuli. In [45] the FH equations were used by Clopton et. a l . to produce tuning curves for single cycles of sinusoidal current stimuli. The model generated 19 tuning curves similar to those obtained from the cochlea of guinea pigs in spite of a number of essential differences between model and experiment. Applications of the FH system involve determining the nature of changes in the model threshold as a function of certain parameters of interest. In the following sections, we present an investigation of the properties of the FH system using techniques of phase-space analysis. Based upon the results of the phase-space analysis, a definition of threshold stimulus amplitude w i l l be proposed, and w i l l be used to produce threshold amplitude l o c i for various current stimuli. The advantages of the phase-based definition of threshold w i l l be shown to be i n terms of increased computational efficiency, and lack of ambiguity. 2.2 NUMERICAL SOLUTION OF THE FH SYSTEM The FH equations were numerically solved by using the predictor-correc-tor algorithm described by Cooley and Dodge [46]. The integration time interval was set at 0.5 |j.s following the recommendations of Butikofer [47] to yield an accurate numerical calculation with a reasonable computation time. A l l calculations were performed using a PDP-11/23 computer. Numerical calculations were checked by comparing data for the time course of V, P^, PR, Pp; m, n, p, h; 1^, IR, 1^, I (Fig. 2.1) in response to a 'standard' stimulus defined by Frankenhaeuser and Huxley [1]. Compari-son with the standard responses i n [1] was qualitative as no numerical data were published. However, the responses appeared to be the same with peak and minimum values occurring at the same times with similar magnitudes as far as could be determined from the graphic data. 20 Fig. 2.1 Standard Response of the FH Model Computed response of the FH model to 'standard* 120 us duration 1.0 mA/cm2 monophaslc current stimulus, a) Transmembrane potential V. b) Ionic current components; sodium, I j ja; potanssium, 1^; u n s p e c i f i c , Ipj leakage, I^. c) Membrane permeabilities; sodium, PN a; potassium, PR; unspecific, Pp. d) Sodium activation variable m, sodium inactivation variable, n; unspecific. activation variable p. 21 2.3 CLASSIFICATION OF THE FH SYSTEM In [48], FitzHugh proposed three classifications for mathematical models of threshold nerve membrane behaviour. The discontinuous threshold phenome-non (DTP) "... in which differential equations with discontinuous functions provide a discontinuity of response as a function of stimulus intensity at threshold ..." [48]. The singular-point threshold phenomenon (STP) was pro-duced by a system of differential equations described by analytic functions with a saddle point present i n the phase space. Such a system has an exact threshold. The third class, the quasi threshold phenomenon (QTP), corre-sponded to a system described by analytic functions with neither a dis-continuous response, nor a sharp threshold. It was shown in [48] and in [49] that the Hodgkin-Huxley model of nerve was of the QTP class. The quasi-threshold nature of the HH equations led to obvious graded responses especially at higher temperatures as illustrated in [50] and [51]. In [51], Cole, Guttman and Bezanilla provided experimental confirmation of the theoretically predicted graded responses with a set of experiments on squid giant axon at temperatures up to 35°C. The Frankenhaeuser-Huxley system i s , like the Hodgkin-Huxley system, of the quasi-threshold variety. The five equations describing the FH system are of the form: dV = F1(m,h,n,p,V) £ - F 2 ( . . V ) . § - F 3 (h ,V) > £ = F 4C„,V), *E - FjCp.V) The set of a l l possible c r i t i c a l points of the five-dimensional system i s 22 given by the use of F± = 0, i = 2,3,4,5 in F^m.h.n.p.V) = 0. The projection of this set of points into the two-dimensional V-m plane in F i g . 2.2 shows that only one intersection (0,0) with ^ = o exists. By the definition of FitzHugh [48], the FH equations thus form a quasi-threshold system. 2.4 PHASE SPACE ANALYSIS OF THE FH SYSTEM The FH system of differential equations expresses the time derivatives of i t s five principal variables, V, m, h, n and p (discussed in detail in a following section) as functions of one another. Solutions to the system would yield expressions for each of the five variables as functions of time, V=<()1(t), m=<t>2(t) ... p=<t>5(t). Such a solution could be interpreted as a point tracing a curve in the space formed by the five variables, with the <t>^ (t) being a parametric representation of this curve. The space formed by the five variables, with time as a paramter, i s known as a phase space. A plane formed by two of the variables, say V and m, would be the Vm phase plane. A curve in the phase space described parametrically by a solution of the FH system is usually referred to as an orbit, path or phase trajectory of the system. 23 A phase trajectory of the syst em w i l l exhibit behaviour which is governed by the locations in the phase space where the time derivatives of the five variables forming the system are zero. Curves in the phase space corresponding to points where the time derivatives of the system variables are zero are known as nullclines or isoclines. Phase trajectories are also highly influenced by points in the phase space where isoclines intersect. These points of intersection, or c r i t i c a l points, may belong to several dif-ferent classes, with each class influencing system st a b i l i t y and excitability in a different manner. For more information on phase spaces, trajectories, and c r i t i c a l points, the reader is referred to [52] by Boyce and DiPrima. 2.4.1 Choice of Phase Space - A Reduced System A reduced system for the Hodgkin-Huxley (HH) equations was proposed by FitzHugh [49] based upon the observation that i n i t i a l l y only two of the four HH variables changed significantly. The remaining two variables were fixed at their resting values and the resulting two-variable system was studied with the aid of the corresponding phase plane. Inspection of responses for the FH variables V, m, n, p, h shows that i n i t i a l l y , as in the HH system, only the V and m variables undergo appreci-able change. The standard responses in Fig. 2.1 i l l u s t r a t e the fact that V (transmembrane potential) and m (sodium current activation variable) are at f i r s t the most active of the five FH variables. Physically, the action potential has a rising phase and a f a l l i n g phase. The i n i t i a l rising phase occurs as the transmembrane potential is driven towards the sodium Nernst potential in response to an increasing permeability of the membrane to sodium 24 V (volts) Figure 2 .2 Quasi-Threshold Configuration of the FH Phase Plane P r o j e c t i o n of set of FH c r i t i c a l points onto the Vm phase plane. Intersection of the locus with dm/dt-0 corresponds to a stable resting point, hence the system i s QTP. 25 ions. The falling phase is due to an increased membrane permeability to potassium ions. This simple physical argument also leads one to the conclu-sion that in i t i a l activity should be confined to the sodium system of any neural model. The m and V variables were thus chosen to form the FH phase plane. With the remaining three variables n, p, h at their resting values, the resulting phase plane was representative of the system for brief periods of time after stimulation. For longer periods of time after stimulation, the same V-m phase plane was used, but with constant values of the variables n, p, h chosen which were more representative of the state of the FH system at the time of interest. 2.4.2 Isoclines and Critical Points in the V-m Phase Plane With the three variables n, p, h held constant at their resting values n , p , h , the loci of points corresponding to F. =0 and F„ = 0 may be o o o 1 2 calculated. From the FH equations with zero stimulus current; or • z (-?NA>GHB<v> - V 0 2 (¥ v>" Vo2(Wv> " 8 * ( V " V ) where GNfl(V) and GR(V) are the Goldman factors for sodium and potassium dV dm respectively. Solving the above equations for — = 0 and — = 0 gives two equations for m in terms of V, dV n 2 -SA(V - V - V Q ^ C V ) - \»Q\W = 0 : m = ^ 1) *w h GM (V) Na o Na 26 dm m „ . dT= ° * m = a (V) + B (V) 2 ) m m Since 0 < m < 1 the negative root in 1) is ignored, and only the positive root i s retained. Since G„ (V) < 0 for V < 0.123 volts, the value of m2 is Na positive. The factor G^a(V) becomes zero for V approximately equal to 123 mV (at 20°C), hence a vertical asymptote of the V isocline i s expected near that value. The denominator of the m isocline term i s never zero, however, i t asymptotically approaches the numerator a (V) as V increases. As V decreas-m es, <xm(V) approaches zero, and thus we expect the m isocline to be a ramp shaped curve starting near zero increasing to some constant value less than 1 for larger values of V. The two isoclines are plotted in Fig. 2.3 as functions of the independent variable V, by regularly increasing V over the interval 0 < v < 120.0 mV and calculating the corresponding m values from 1) and 2). With n, p, h held constant at their i n i t i a l values, there are three intersections of the isoclines, and hence three c r i t i c a l points in the V-m phase plane. The point (0.0 mV, 0.0005) is a stable improper node, the point (17.0 mV, 0.055) is an saddle point, and the point (120.0 mV, 0.98) is another stable improper node. The values of V and m corresponding to the c r i t i c a l points were deter-mined by l i n e a r i n t e r p o l a t i o n to locate the points where m = V. To more accurately determine the location of the stable resting point, n, p, h were held constant at their resting values and the response to a null stimulus was calculated. The resulting values of m and V quickly approached those corres-27 Figure 2.3 Resting State Isoclines of the Vm Pahse Plane Resting state isoclines dV/dt-0, dm/dt-0, of the reduced Vm phase plane. 28 ponding to the stable resting point - these values were (-0.16 mV, 0.00048). Similarly, for the stable point near 120 mV, the equations were solved with n, p, h constant at their resting values, starting with V=120 mV, m=0.98. The computed response converged quickly to (119.3 mV, 0.994). For the saddle point near 17 mV, the sign of was altered to transform the point into a stable critical point and the computed response from an in i t i a l value in the neighbourhood of the point was observed to converge to (17.0 mV, 0.0544). 2.4.3 Threshold Characteristics of the FH System Although the complete FH system is quasi-threshold in performance, for certain restricted parameter ranges it exhibits singular-point threshold behaviour. The critical parameter ranges are illustrated in Figs. 2.5-2.8 and are discussed below. The saddle point at (17.0 mV, 0.0544) gives the FH system its sharp threshold characteristics for certain parameter ranges. With one exception, a l l trajectories eventually diverge from the neighbourhood of the unstable saddle point. The single convergent trajectory separates two families of divergent trajectories of which one group will veer towards the stable resting node while the other group will head towards the stable excited node. These two distinct trajectory groups represent supra-threshold and sub-threshold trajectories of the FH system. The single convergent trajectory, or separatrix, may be calculated by numerical reverse integration from a starting value in the neighbourhood of the saddle point. The separatrix resulting from such a calculation is shown in Fig. 2.4 along with a 29 Figure 2.4 Keat ing State Phase Plane D e t a i l Res t ing state phase plane near the unstable saddle p o i n t . The separatr ix (dotted) d i v i d e s sub- threshold ( lef tmost ) and supra- threshold (rightmost) phase plane t r a j e c t o r i e s . 30 representative trajectory from each of the supra-threshold and sub-threshold groups. The presence of this saddle point, and hence the STP behaviour of the FH system, is transient for the full set of variables and for any reduced system which includes the effects of recovery from the excited state. 2.4.4 Phase Plane Variations with n, p, h Although the presence of a stable excited node is indicated by this reduced phase analysis, calculated solutions of the full FH equations show no such stable excited response. The excited state of the f u l l system is transient and the computed transmembrane potential returns to its original stable resting value. The reduced phase plane, however, represents only that part of the five-dimensional phase space formed by the intersection of the space with the surface defined by the i n i t i a l values of n, p, h. Hence variations in the other three variables n, p, h must be considered to under-stand the transient nature of the excited state. Assuming the FH standard response to be fairly typical, the following variations in the parameters n, p, h during the course of an action potential were noted: 0.0269 < n < 0.56 0.0049 < p < 0.188 0.02 < h < 0.8249 To observe the effects of changing n, p, h upon the phase plane, isoclines were re-calculated for values of each variable within the above ranges. 31 2.4.4.1 Potassium Activation Variable n As n increases from its in i t i a l value of 0.0268, the magnitude of the potassium current component also increases. In Fig. 2.5 are a series of level curves of dV/dt = 0 for various values of n. The m isocline is not altered by changes in n. The general effect of increasing n is to raise the activation threshold through a translocation of the unstable saddle point to a position 'higher' on the m isocline, and to 'lower' the stable excited node. 2.4.4.2 Unspecific Current Activation Variable p As p increases from its i n i t i a l value of 0.0049 to its maximum near 0.188, the unspecific current component increases. A series of plots of the V isocline for various p values are shown in Fig. 2.6. Increasing p has the effect of causing a coalescence of the stable resting node and the unstable saddle point until, i f p is large enough, the two critical points vanish. In this case only the stable excited node remains, and every trajectory becomes supra-threshold. The saddle point vanishes for p > 0.13 (approximately). 2.4.4.3 Sodium Inactivation Variable h The resting value of h is the maximum value attained by h during the course of an action potential. Decreasing h from its resting value decreases the magnitude of the sodium current and, i f h is small enough, turns the current off. In Fig. 2.7 a series of V isoclines is plotted in the phase plane for various h values. 32 Figure 2.5 Affects of Changing 'n' on the Vm Phase Plane The affect of changing potassium activation variable, n, on the phase pi The m isocline in unaffected. 33 Figure 2.6 Affects of Changing 'p' on the Vm Phase Plane The affect of changing unspecific activation variable, p, on the phase plane. The m isocline i s unaffected. 34 Figure 2.7 A f f e c t s of Changing * h ' on the Vm Phase Plane The a f f e c t of changing sodium a c t i v a t i o n v a r i a b l e , h , on the phase p i The m i s o l c i n e i s u n a f f e c t e d . 35 As h decreases, no change in the position of the resting node occurs. However, the excited node i s shifted to the le f t along the m isocline and the saddle point is shifted to the right along the m isocline. These two nodes ultimately coalesce at some value of h near h = 0.02, then both nodes vanish for smaller h values leaving the resting node as the only c r i t i c a l point in the system. 2.4. A.A Sodium and Potassium Permeability From Fig. lc) and d) values of the variables n, p and h can be determined at the instants corresponding to peak membrane permeability to sodium and potassium. The phase plane i s shown in Fig. 2.8 as i t appears for maximum sodium and potassium permeability. At peak P^a the phase plane changes are not large. A slight increase in threshold i s observed due to the small change in saddle point location. At dm the same time the stable excited node moves le f t along = 0. At peak P^ the phase plane i s d r a s t i c a l l y a l t e r e d . Both the stable excited node and the saddle point have vanished, leaving the stable resting node as the only c r i t i c a l point i n the phase plane. At peak P^ the system i s thus completely un-excitable in a state corresponding to the absolute refractory period of the nerve, and w i l l remain so until the V isocline again intersects the m isocline to produce a stable excited node. When the V dm i s o c l i n e does re-intersect - 3 — = 0, the threshold w i l l i n i t i a l l y be quite d t high, since the saddle point w i l l be located high on the m isocline, but w i l l decrease as the saddle point approaches i t s resting location i n the phase-plane - this corresponds to the nerve's relative refractory period. 36 2.5 DETERMINING THRESHOLD STIMULUS PARAMETERS Given the sharp threshold characteristics of the FH system when a saddle point i s present in the V-m phase planes, i t is of interest to be able to determine the minimum stimulus amplitude which w i l l e l i c i t an action potential. In [53] Frankenhaeuser computed a series of responses to short duration monophasic current stimuli "... to find a threshold strength for regenerative a c t i v i t y " . To do th i s , a stimulus amplitude was changed in successive runs such that i t converged to the limiting threshold strength. Apparently Frankenhaeuser identified stimuli as supra-threshold by observing the com-puted action potential - i f the computed potentials did not show an active response (by demonstrating the typical 120 mV action potential) they were classified as sub-threshold. The intuitive meaning of threshold stimulus amplitude i s clear, but more precise definitions are required to automate threshold detection and to remove descriptive ambiguities which can arise when applying the intuitive threshold amplitude definition in [51]. FitzHugh provided two analytical definitions of threshold stimulus amplitude, the level definition and the inflection definition. While these definitions were formulated for the Hodgkin-Huxley system, they are valid for the Frankenhaeuser-Huxley equations as well. The l e v e l definition of threshold requires a c r i t i c a l value, V , of the transmembrane potential V to be selected, and to be used in conjunction with a stimulus-response curve to define threshold stimulus amplitude. The stimulus-response curve i s a continuous curve of peak transmembrane poten-37 Figure 2.8 Vm Phase Plane at Peak PN a and Peak PR The phase plane at the instant of maximum sodium permeability, PN_ (h - 0.76, p - 0 .03, n - 0 .13) , and maximum potassium permeability, PR(h - 0.025), n - 0.555, p - 0.19). 38 t i a l , V , as a function of stimulus amplitude for some particular stimulus waveform. A set of stimulus-response curves for FH standard monophasic square stimuli at two temperatures are shown in Fig. 2.9. According to the level definition, the threshold stimulus amplitude i s the stimulus amplitude at which the stimulus response curve crosses the c r i t i c a l l e v e l V . c Arbitrarily choosing Vc=60 mV, this gives a threshold of about 0.75 mA/cm2 at 20°C, and 0.84 mA/cm2 at 15°C for the FH standard 120 us monophasic pulse. The value of threshold amplitude obtained in this manner i s only slightly affected by the choice of Vc as long as the maximum slope of the stimulus-response curve i s large. For a DTP or STP system, the slope w i l l always be large. For a OTP system the responses may be graded, and hence the slope values may be more moderate. The inflection definition of threshold requires that one calculate the slope of the stimulus-response curve at various stimulus amplitudes. The stimulus amplitude corresponding to maximum stimulus-response slope i s defined to be the threshold amplitude. Although requiring more computation than the level definition, the inflection definition has the advantage of giving a unique threshold amplitude even for stimulus-response curves showing a graded, rather than sharp, threshold response. For a given stimulus waveform, to determine the threshold amplitude, a series of computations are required to obtain a number of points on the stimulus response curve. If the level definition i s being used, computation w i l l cease when the peak response exceeds the previously established V le v e l . If one i s using the inflection definition of threshold then the 39 Figure 2.9 Stimulus-Response Curves f o r 120 ws Monophasic Square S t i m u l i at 2 0 ° C , and 1 5 ° C . Threshold amplitudes determined by the l e v e l d e f i n i t i o n of threshold with V £ ^ 0 mV; 0.75 mA/cm2 @ 2 0 ° C , and 0.84 mA/cm 2 @ 1 5 ° C . 40 threshold amplitude has been reached when the next stimulus amplitude results in a lower stimulus-response slope. The inflection definition was used by FitzHugh [50] to calculate thresh-old amplitudes for anodal break excitation in the Hodgkin-Huxley system. A l e v e l d e f i n i t i o n was used by Butikofer and Lawrence [44] with V c =80 mV to produce threshold charge curves for biphasic square stimuli. While both definitions were adequate in each of the applications, the definitions are computationally inefficient and, in extreme cases, ambiguous. Ambiguities can arise when trying to determine thresholds for short duration biphasic pulses. Such pulses can generate transient membrane potentials of almost arbitrary magnitude while e l i c i t i n g an overall sub-threshold response and subsequently confusing the stimulus-response picture. An example of such a response i s shown in Fig. 2 . 1 0 . To avoid these ambigui-ties for short biphasic stimuli i t would be necessary to exclude from consideration responses during the stimulus, and to note the maximum response upon conclusion of the stimulus. Computational inefficiencies result due to the lack of a positive definition of sub-threshold response. Thus the threshold program may continue to compute a stimulus response when, i n fact, there i s no possibil-i t y of that response generating an action potential. One normally would select a time limit and state that a response which was s t i l l not active at the end of the time limit was sub-threshold. Shortening the allowed time interval to improve computation effiency, however, w i l l have the effect of a r t i f i c i a l l y inflating the threshold amplitude. 41 8 0) > 0) o Q . c 01 fg c 8 o'.io Tiae oTaJ oTao (•i Ill-seconds) 'o.oo 0 . 4 0 Figure 2.10 Sub-theshold Response to 5 ps Biphasic Square Stimulus 0.0 ps Phase Delay Subthreshold response to a 34 mA/cm2 5 us wide biphasic square stimulus with 0.0 us phase delay. Note that threshold, according to the level definition, i s exceeded in this example. 42 2.5.1 The Separatrix Definition of Threshold The phase plane and separatrix of Fig. 2.4 suggest an unambiguous and computationally efficient threshold definition which may be applied as long as the saddle point and associated separatrix exist in the V-m phase plane. I f , at the completion of a stimulus of duration t the resulting phase trajectory has not escaped to the supra-threshold region of the phase plane as defined by the separatrix, then the stimulus i s sub-threshold. Thus the problem of threshold definition becomes one of determining on which side of the separatrix a point l i e s which corresponds to the phase plane value of the trajectory at stimulus end. In the case of an exponentially decaying stimulus, some convenient definition of stimulus end must be employed. To implement this threshold definition, a least-squares surface was fitted to a series of separatrices corresponding to V-m separatrix curves for various values of the sodium inactivation variable, h. This surface gave the sodium activation variable m as a function of V and h. A stimulus of dura-tion t w i l l give r i s e to a point (Vt, mt) in the phase plane with a corre-sponding h value of h^, when the stimulus has just ended. If the separatrix i s given by m =m (V,h), then a stimulus i s defined as being supra-threshold s s i f and only i f mt > VVV Hence, at stimulus end, the nature of the response can immediately be deter-mined. This definition of threshold is insensitive to large transient fluctuations in V since i t is concerned only with the value of V, m, and h immediately after the stimulus i s over. It offers a significant improvement 43 in computational efficiency since i t i s unnecessary to perform a protracted computation to ensure that no delayed activity i s possible. The separatrix definition does not provide a definition of the threshold amplitude as do the level or inflection definitions. However i t provides a means of determining threshold amplitude to an arbitrary degree of accuracy. Let e>0 be the accuracy parameter; i f , lmt " m s ( V t , h t ) l < £ a n d m t > W V then the stimulus amplitude giving r i s e to (Vt» m^ ) i s the threshold amplitude. If i - i (V^.li ) < 0 then the stimulus i s sub-threshold, i f t s t t nr - m (v^,h ) > e , then the stimulus i s supra-threshold. In practice, e t s t ' t ' would be bounded below by some f i n i t e non-zero value which would depend upon the degree of influence of the neglected n and p. This treatment of theshold separatrix as a series of level curves in the V-m plane, parameterized by the inactivation variable h, ignores the effect of the variables n and p upon threshold and determines threshold by using the f u l l sodium system. With the exception of long duration stimuli, or repetitive stimuli, neglecting the influence of n and p upon the separatrix is an acceptable simplification. Threshold amplitude l o c i for some typical biphasic and monophasic square pulses determined using the separatrix definition are shown in Fig. 2.11. Although only results up to 100 (is duration are reported, i t was observed that threshold amplitudes of pulses up to 500 us i n duration could be determined to 1% accuracy. Biphasic pulses longer than 1 ms gave consistently inaccurate results, that I s , amplitudes reported as supra-44 Figure 2 .11 Threshold Amplitudes of Biphasic Stimuli Determined by the Separatrix Definition Threshold amplitude determined by separatrix definition as a function of biphasic current pulse width with phase separation as parameter. Infinite separation corresponds to monophasic pulse. 45 threshold actually generated no action potentials, however i t was possible to determine to within 1% the threshold amplitude of a monophasic square pulse 1 ms in duration. The sources of error in this method of threshold amplitude calculation include; neglecting n and p effects on threshold, using an approximate surface to f i t the separatrix, and neglecting separatrix changes following the stimulus. It does not seem possible to directly calculate the error introduced into the threshold amplitude, but i t can be inferred that the error is less than 1%. A l l threshold amplitude values obtained by the separatrix definition were tested by computing the response of a f u l l FH system to the particular stimulus with the threshold amplitude. A l l such threshold stimuli were observed to generate action potentials and i f the amplitude was reduced by 1% - a sub-threshold response was noted. In Fig. 2.12 the ratio of the time taken to compute threshold using the level definition to the time taken using the separatrix definition i s plotted for the case of a monophasic square stimulus. The greatest improvements are observed at lower pulse durations, the separatrix method being over 60x faster for a 5 us wide pulse. At 100 us pulse durations the improvement factor i s a more modest 5x. Simulation of 1 ms of membrane activity requires approximately 120 seconds on a PDP 11/23 with an extended instruction set for a model coded in PASCAL. The parameters governing the numerical integration algorithm are set to yield highly accurate results, the integration interval being 0.5 us with a corrector tolerance of 10-lf. Increasing either of these parameters would decrease the execution time at the expense of numerical accuracy. 46 8 Figure 2.12 Comparison of Threshold Calculation Times for Separatrix and Level Definitions Ratio of time taken to calculate threshold amplitude of monophasic pulse using level definition, to time taken using separatrix definition, as a function of pulse width. 47 In calculating threshold stimulus amplitudes, execution time is influenced by the i n i t i a l l y provided supra-threshold amplitude, and, as previously mentioned, by the simulation time limit imposed during each cycle of the threshold iteration process. When using a level definition of threshold at 20°C, one must generally choose a limiting time no less than twice that of the stimulus duration unless the stimulus i s quite brief in which case a limiting time more than twice the stimulus duration may be required. Temperature w i l l influence the limiting time for any given stimulus• Specifically, for a monophasic current pulse of 100 us duration, starting with threshold bracketing amplitudes of 1.0 mA/cm2 and 0.0 cm2, with a limiting time of 750 us per iteration the following execution times resulted: - separatrix definition: 93 sec. - level defintion: 589 sec. Both definitions gave a threshold amplitude of 0.8672 mA/cm2. 2.6 DISCUSSION A measure of the "sharpness" of the Frankenhaeuser-Huxley QTP adapted from that of FitzHugh [48] is given by: |V(I2;t2) - V ( I1; t2) | |V(I2;t1) - V d ^ t ^ l where I2 and I± are current pulses of duration t j , respectively producing 48 " a l l " and "none" responses, V(I;t) is the computed response to I at time t , and t2 is the time of maximum response to stimulus I j - The larger this ra t i o , the more sharp is the QTP. In [48] FitzHugh reported a sharpness of 9.8*105 for the Hodgkin-Huxley system based upon unpublished data of Cole. Similar calculations for the FH system show that for monophasic Dirac current pulses, slightly supra-thres-hold and slightly sub-threshold pulses yield transmembrane potentials of 29.7871 mV and 29.7872 mV respectively - a difference of 10 nV. The differ-ence in V for the two cases measured at the time corresponding to the peak response of the supra-threshold I2 i s 103.3 mV. This establishes a lower bound of 1*107 for the sharpness measure of the FH system. A similar sharp-ness calculation for 100 us current pulses yields a value of 1.1*107. Since the FH system exhibits a sharpness parameter two orders of magni-tude greater than that of the HH system, one would expect to see less res-ponse gradation by the FH system. This i s , in fact, the case as illustrated in F ig. 2.13. This figure shows stimulus response curves of the FH system for a 120 us monophasic pulse at temperatures of 2.5°C, 20°C and 40°C. It must be noted that a temperature of 40°C i s a significant extrapolation beyond the original 2.5°C to 22°C data base for which the various Q^Q'S W E R E measured. The rounding of the 40°C curve is due to the fact that the Q10 parameters for the sodium inactivation system are larger than those of the sodium activation system thus increasing the rate of accommodation at higher temperatures. The curve exhibits no gradation, however, in comparison with the HH system at similar temperatures [50], [51]. It i s interesting to note that Increasing temperature past 40°C eventually results in a transformation 49 Figure 2.13 Stimulus-Response Curves for 120 us Monophasic Square Stimuli at 2.5° C, 20° C, and 40°C The curves exhibit l i t t l e response gradation even at high temperatures. This i s a reflection of the large 'sharpness' parameter associated with the FH model. 50 of the FH system to a singular-point phenomenon. This, however, occurs at a temperature of 375°C and hence is of no practical value. 2.7 SUMMARY The Frankenhaeuser-Huxley model of myelinated nerve has been shown to be a system of the quasi-threshold type for physically meaningful temperatures, with a sharpness parameter two orders of magnitude greater than that of the Hodgkin-Huxley system. The sharpness of this system is such that graded stimulus-response curves are not observed, even at relatively high temperatures. The excitability of the FH system, investigated using phase plane analy-s i s , was shown to depend upon the location of a saddle point and i t s associated separatrix in a reduced phase plane. By f i t t i n g a least-squares surface to the separatrix as a function of the i n i t i a l l y most active FH variables V, m and h, i t was possible to automate a separatrix threshold definition which resulted in a significant reduction of computation time when iterating to a value of threshold stimulus amplitude. This separatrix definition of threshold i s unambiguous in cases where the level and inflec-tion threshold definitions f a i l . Although not directly applicable to the determination of experimental thresholds this definition provides an algor-ithm that can be used to extract computer model thresholds for comparison with experimental data. 51 3. OPTIMIZATION OF NEURAL STIMULI  BASED UPON A VARIABLE THRESHOLD POTENTIAL 3.1 INTRODUCTION The problem of determining the minimum stimulus amplitude required to e l i c i t a nerve action potential - the threshold problem - is one of consider-able practical significance. With knowledge of the l o c i of minimum ampli-tudes as functions of various stimulus parameters and waveform types, one can establish a heirarchy of stimuli ordered by relative comparisons of threshold amplitudes. We refer here to 'amplitude' as the current parameter which i s a constant in the mathematical expression of time-varying stimulus current. Although stimulus current may also be defined in terms of peak or average current, energy or charge - we propose to discuss threshold with respect to a parameter which may be controlled directly at the simulator, rather than discussing threshold with respect to parameters which would be calculated from, or functions of, the stimulus amplitude setting. For any given stimu-lus, the threshold locus may be inspected to select a range of parameters minimizing stimulus amplitude or charge per phase of the stimulus while s t i l l ensuring generation of a nerve action potential. Minimizing stimulus ampli-tude minimizes an amplitude dependent tissue damage component due to ohmic heating [11], while minimizing stimulus charge maximizes the l i f e of the stimulator's power source and minimizes a charge proportional damage com-ponent [55], [9]. It i s variously reported in the literature that the tem-perature increases responsible for heating damage are large enough to cause local blackening of the skin [12], and small enough to be negligible [56]. 52 Such threshold amplitude or charge l o c i are, however, analytically elusive. The familiar strength duration (SD) curve provides a useable threshold amplitude function for the case of monophasic square constant-current stimuli, but few other deterministic expressions of threshold ampli-tude have been previously reported. The emphasis of this paper w i l l be on developing a simple analytic expression of threshold amplitude for the case of monophasic and biphasic constant current stimuli. The l o c i of threshold amplitude (as a function of stimulus duration at some constant phase separation) produced by the afore-mentioned function, w i l l be compared to numerical results from the Frankenhaeuser-Huxley (FH) model of myelinated nerve and to the threshold amplitude l o c i of the classical SD curve. The expression for threshold cur-rent amplitude w i l l also be used to describe the stimulus waveform which minimizes tissue damage due to both thermal power dissipation and charge injection. In [16], McNeal presented a model emulating the response of an eleven-node myelinated fibre to constant current stimuli. The intention of McNeal's work was to develop a modelling system where the effects on neural threshold of arbitrary electrode configurations and fibre geometries could be assessed. With emphasis on the geometric aspects of neural threshold, McNeal used passive Frankenhaeuser-Huxley systems at each of the eleven nodes. The membrane impedances were fixed at values corresponding to the quiescent state of the FH system and were not allowed to vary with transmembrane potential. The emphasis in this paper i s away from the geometric influences of electrodes and fibres, so well covered by McNeal, and more towards the 53 implications that membrane non-linearities have with regard to threshold and strength-duration relations. Hence, our approach i s to calculate responses to current stimuli of a single fu l l y active FH node which has a l l physical parameters normalized with respect to unit surface area, ignoring the effects of electrode geometry. As pointed out by McNeal in [16], this amounts to an assumption of a stimulus applied to a single node which i s valid in the case of stimulation with a microelectrode. Experimental results in [44] demon-strate that single node FH model data agree closely with threshold data from transcutaneously applied stimuli, indicating that the assumption may be extendable to stimulus scenarios other than that of percutaneous micro-electrodes. 3.2 THE CLASSICAL STRENGTH-DURATION CURVE The strength-duration curve, an expression of threshold amplitude for constant current monophasic square waves, f i r s t derived by Lapicque in 1907 [57], takes the form WT> • w * 1 - e " T / T > ( 1 ) where T i s the pulse duration, T the membrane time constant and I ^  the rheobase current. As pointed out by Noble and Stein [58], although the constants in this equation may be functions of the stimulation scenario (electrode geometry, distance from the nerve etc.), the form of the relation is independent of stimulation scenario thus making this relation widely applicable. 54 The derivation of this equation i s predicated upon two assumptions: 1) the nerve response may be approximated by a passive RC network, and 2) the nerve threshold may be taken to be some constant AVc above resting potential. The assumption of passive RC response i s recognized to be invalid i f the transmembrane potential becomes too large (> 20 mV), or i f accommodative effects become significant as in the case of long duration stimuli. The constant threshold approximation was proposed by Hodgkin and Rushton in 1946 [59], based upon a core conducting cable model. Given a constant current stimulus of amplitude I, assuming an RC mem-brane, the membrane potential above resting potential becomes V(t) = I.R.(1 - e"t / T) (2) where R is the membrane resistance, C the membrane capacity and T = RC the membrane time constant [60].- For a constant current monophasic pulse of duration T, assuming a constant threshold potential of AVC> the threshold amplitude, It n» i s given by It h = AVc/{R.(l - e"T / T)} (3) The rheobasic current, I , , i s defined as lim I , (T), and thus I , = AV /R ' r h ' „ thv " rh c with the resulting SD curve given by \* - W*1 - rI/T> <4> Since t = RC may be calculated from known membrane constants for R and C, only I remains to be determined. This would usually be done directly by noting the threshold stimulus amplitude approached by long duration pulses 55 rather than by measurement of the constant potential threshold AV^ and subse-quent calculation of the rheobasic current from I , = AV /R. rh c Using a threshold amplitude locus numerically obtained from the FH model as a basis for comparison (at 20°C, with 1/R = 0.0303 m h o » c m2, C = 2.0 uF/cm), i t can be seen that the hyperbolic strength-duration relation does not agree very well with the corresponding FH threshold amplitudes at low pulse widths (Fig. 3.1). Agreement i s good for pulse widths large with respect to x, since the asymptotic I ^ value was chosen to match the observed threshold current for long duration pulses. 3.2.1 Limitations of the RC Strength-Duration Curve Apart from increasing deviation from the FH model results with decreasing pulse widths and the i n a b i l i t y to predict long duration thresholds due to the exclusion of accommodation influences, the RC strength-duration curve cannot be applied in the case of biphasic stimuli. This is because the constant threshold potential of the passive RC network modelling the nerve may be exceeded only during the leading phase of a biphasic stimulus. The biphasic thresholds so determined w i l l be independent of phase separation, a situation which is known not to be true [4], [6] and [15]. In the case of long duration stimuli, the accommodation insensitive RC network i s not a good neural model. In the case of monophasic stimuli and biphasic stimuli with phase separation greater than about 20 us, the trans-membrane potential does not greatly deviate from the form produced by a passive RC network. This can be determined by observing responses of the FH model to the stimuli in question. The source of the strength-duration inade-56 Figure 3.1 Model Strength-Duration Curves Compared to FH Model Threshold Amplitudes 57 quacies for short pulses and biphasic pulses i s thus the assumption of con-stant potential threshold. In [15], i t was shown that a threshold separatrix exists i n the phase space of the FH system. This separatrix, which divides sub-threshold phase trajectories from supra-threshold trajectories was shown to be a strong func-tion of the sodium system of the FH model. In particular, Fig. 2.4 (Fig. 4 of [15] demonstrates that the threshold separatrix i s not a constant function of V, thus implying that the threshold potential function i s not constant either. If the assumption of constant threshold potential was correct, then the separatrix of Fig. 2.4 would appear as a plane in FH phase space, inter-secting the V-m phase plane as a vertical line independent of m, passing through the FH unstable saddle point• 3.3 THE THRESHOLD POTENTIAL FUNCTION If the assumption of constant potential threshold i s not adequate for certain pulse width ranges and for certain types of stimulus waveforms, then what form would a variable potential threshold take? The FH model of myeli-nated nerve was used to provide some insight into this question for mono-phasic and biphasic stimuli. Calculations with the model reveal that there i s a value of transmem-brane potential which must be attained by stimulus end in order to e l i c i t an action potential. However, the potential which must be exceeded is not a constant but instead appears to be a function of pulse width with short dura-tion pulses requiring higher potentials than longer stimuli. 58 3.3.1 Monophasic Stimuli The locus of these FH system threshold potentials for monophasic stimuli of various durations i s shown in Fig. 3.2. These values were obtained by calculating threshold stimulus amplitudes and then recording the transmem-brane potential attained by the FH system at the end of the stimulus. At any given stimulus duration T, the threshold stimulus amplitude was determined such that a 0.1% reduction in stimulus amplitude produced a sub-threshold response. This locus of V ^ i s clearly a function of stimulus duration T, and i s not constant. As T increases V ^  approaches some asymptotic value, and as T->0, V ^  takes on some other f i n i t e value. These trends suggest a functional relation of the form Vt h(T) = a.e"bT + c (5) The constant 'c' corresponds to the constant threshold potential of the classical SD curve, 'a' is related to the displacement of V required to bring the i n i t i a l phase point of the FH system to rest on the threshold separatrix [15], while 'b' determines the rate of decay from the maximum threshold potential at T=0, to the asymptotic value 'c'. 3.3.2 Biphasic Stimuli Also shown in Fig. 3.2 are the threshold potential l o c i for biphasic square stimuli with various phase separations. For consistency of definition with the monophasic potential thresholds, the biphasic potential thresholds were defined to be the transmembrane potentials attained by the FH system at the end of the leading phase of the stimulus. This i s consistent with 59 Figure 3.2 Threshold Potentials for Biphasic Square Stimuli from the FH Model Threshold potential increments (above resting potential) required for excita-tion at various pulse widths, for monophasic and biphasic square constant current stimuli. 60 consideration of the monophasic pulse as a biphasic stimulus with the t r a i l -ing pulse i n f i n i t e l y removed. The biphasic threshold potential l o c i for pulses with phase separation no less than about 20 us, can be represented by functions of the same form as the monophasic case Vt h(T) = a.e"bT + c (6) The parameters a, b, and c are constant for any given pulse separation and are numerically determined by curve f i t t i n g . In general, a l l three para-meters increase as phase separation decreases. 3.4 THRESHOLD CURRENT AND CHARGE For monophasic stimuli and biphasic stimuli with phase separations not less than about 20 us, a potential threshold function requiring three constants with one exponential term may be used. The resulting expression for threshold current as a function of pulse width at a particular phase separation i s V e Ir h *( 1 " e 5 t h (1 - e T / t) (1 - e T / t) using the variable threshold in conjunction with the assumption of an RC membrane. The constant 'b' is the exponent which appears i n the V ^  expres-s i o n , while x i s the membrane time constant. The constant I , is identical rn to the previously defined rheobasic current of the classical strength-duration curve and is determined from lim I , (T) with _ th I r h - c/R (8) 61 where ' c' is the asymptotic potential threshold constant. The parameter 1^, the theta current, i s related to the threshold potential of a Dirac delta function stimulus, and is determined from constants 'a' and 'c' in the expression for exponential potential threshold. For an instantaneous thresh-old current pulse, some charge Q is injected with a resulting transmembrane potential of AVQ = QQ/C. The I. is defined to be o o o I0 = AVQ/R (9) The threshold current predicted from this form of strength-duration curve is shown in Fig. 3.3 along with the threshold current values determined from the FH model. The corresponding threshold charge function can be determined from Qth = It h(T)*T (10) These calculated threshold charge values w i l l over estimate the actual FH charges due to the 1^(1) results being greater than the corresponding FH threshold amplitudes. As can be seen from Fig. 3.3, the RC model in conjunction with the exponential threshold potential yields reliable results only up to pulses of about 100 us in duration, with phase separations greater than 20 us. 3.5 STIMULATION WITH MINIMUM POWER One of the most interesting results derived from the constant threshold SD relation i s the minimum power waveform of Offner [61]. Using a passive RC network with a constant potential threshold, Offner showed that for a mono-8 *1 62 8 9-=8 r 8 loo oTio oTao oTso . o.40 Pulse Width (•illi-seconds) O.HO 8 8 r 8 rr oTii o'.eo o'.ao o'.40 Pulse Width (•!Ill-seconds) •bToo o.to Figure 3.3 RC-Varlable Threshold Potential Strength-Duration Curves Compared to FH Model Threshold Amplitudes Threshold current amplitudes calculated by the FH model compared to the threshold amplitudes determined from a passive RC model using in conjunction with a variable threshold potential function. 63 phasic stimulus the waveform of duration T minimizing generated heat while s t i l l exciting the nerve has the form of a rising exponential current 2 ( A V t/x ^ f n e r ^ - T/x - T / T , R(e - e ) with x - RC the membrane time constant, and AV the constant threshold ' c potential. Assuming a variable potential threshold and a damage component propor-tional to power, as did Offner, one can minimize T 2 H = k / I (t)dt k = constant, (12) o for a stimulus of duration T, subject to the constraint that the threshold condition i s met within time T, while the transmembrane potential evolves according to -tlx t . V(t) - e-g— • / I(s)es / Tds (13) o This i s accomplished by introducing a Lagrange multiplier, X , to minimize /T( k l2( t ) - £ e- T / T et / TI ( t ) ) d t (14) o subject to the constraint that V(T) = Vt h(T) = a » e "b T+ c (15) From the Euler-Lagrange equation of variational calculus, 64 with the multiplier \ being eliminated to yield R(e - e ) I . (t) = I- e"b T« csch(T/x) • et / T + I... (T) (18) mmv ' R \ ' / Offnerv ' The expression for heat generated by a stimulus is minimized for pulses of i n f i n i t e duration. This, however, is beyond the recognized limits of the simple RC model, and hence any value of I . obtained by evaluating v ' ' min J & lim I . (T) should be cautiously interpreted. T+- m i n The assumption of a potential threshold with three constants yields a minimum power waveform which has the same rising exponential form of the Offner minimum, but with a larger i n i t i a l value. As pulse width T grows, the variable threshold minimum approaches ^offner s*n c e ^-ar8e ^ corresponds to a range of pulse width where the constant potential threshold assumption is more applicable (until the effects of accommodation become significant). The largest differences in the two waveforms occur for small T values. As TO the i n i t i a l I . value exceeds that of I,^,-,- by a factor of (1 + a/AV ) . min Offner J c The constant 'a', parameterized by phase separation, ranges from 9.3 mV for the monophasic case, to 17.2 mV for a biphasic pulse with a 20 us phase separation, corresponding to i n i t i a l current amplitude increases by factors of 1.7 and 2.3 respectively. The variable threshold potential l o c i used to derive the power mini-mizing Offner-like stimuli were produced with reference to biphasic square 65 waveforms and hence may not be applicable to the case of exponentially rising stimuli. Calculation of threshold potential values for the rising exponential waveforms reveals that for short duration stimuli the threshold potentials are virtually identical, while at longer stimulus durations, slight differences in threshold potentials were observed. Quantitatively the following differences were noted: - for 5us duration monophasic stimuli: square stimulus - V , = 29.45 mV i . „ _ n i n th } AV = 0.01 mV exponential stimulus - V ^  = 29.44 mV - for 500 us duration monophasic stimuli: square stimulus - V , = 20.68 mV i , _ th } AV = 1.2 mV exponential stimulus - V ^ = 21.85 mV - for 5 us duration 20 us phase separation biphasic stimuli: square stimulus - V , = 38.74 mV i . „ n T . n th f AV = 0.02 mV exponential stimulus - V n = 38.72 mV - for 500 us duration 20 us phase separation biphasic stimuli: square stimulus - V, = 22.43 mV l »TT -i -i-? ™ M th } AV = 1.17 mV exponential stimulus - V ^  = 23.60 mV At low stimulus durations for a l l phase separations, the square wave threshold potentials slightly exceeded those of the exponential stimuli, while at long stimulus durations for a l l phase separations the square wave threshold potentials were less than those of the corresponding exponential stimuli. This result i s as expected since at short durations the exponential 'spike' causes the FH system to have a larger dV/dt at stimulus end than does a square wave; while at long stimulus widths the i n i t i a l low amplitude of the 66 exponential stimuli means that the sodium activation variable rate of increase, dm/dt, i s low thus implying that a larger threshold potential is required to excite the FH system [15]. The slight differences between threshold potential l o c i for the two types of stimuli should not be significant enough to cause concern over use of the square stimulus potential thresholds to derive the form of the exponential waves. In addition to the stimulation damage caused by heat dissipation in the tissue, there i s a possible charge proportional damage component [9], which in the case of a charge balanced biphasic stimulus would be in terms of charge per phase of the stimulus. Let H be the heat related damage component, Q be the charge related damage component and y be some dimensionless parameter between 0 and 1. Define D(T) to be a function of pulse duration parameterized by phase separation which expresses damage due to stimulation in some normalized dimensionless units as 3.6 STIMULATION WITH MINIMUM CHARGE D(T) - YH + (1 - Y)Q (19) with H / I2(t)dt (20) o T and Q / I(t)dt (21) o 67 In this damage function the value of y would be set to reflect the relative degree of importance placed upon the two damage mechanisms. Minimizing this new function subject to the same threshold constraint as before gives V (T) W T > - ( - ^ - + 1 C l - e-T^))csCh(T/,) • e ^ + \ (22) where vt n( I ) Is the previously defined threshold potential function. Substi-tution of this stimulus waveform into the damage function gives composite damage as D(T) = YxA2(T)csch(T /T) + ( 1 (23) with V (T) A(T) = ( - ^ — eT / 2 T + (LZJL) sinh (T/2x) ). (24) The minimum of this expression i s expected to occur at some f i n i t e non-zero value of T, parameterized by y and phase separation, since minimum charge occurs as T+0 [3], and minimum heating occurs as T-><» [61]. Calculation reveals that the minimum of the damage function occurs when 2T = A(T)coth(T/x) (25) The resulting expression i s awkward and i s best solved numerically. With y selected as 0.5 (equal weight to heat and charge damage components) the following damage minima were obtained for given phase separations: - for a monophasic square pulse, minimum damage at a pulse width of about 400 us 68 - for a biphasic square pulse, with 20 as phase separation, minimum damage at about 500 us pulse width. The locus of damage minima as a function of y is shown in Fi g . 3.4 for both of the above pulses. Since the relation of the parameters determining vtj1(T) in terms of phase separation i s not e x p l i c i t , i t is only possible to minimize D(T) on a constant separation isosurface, and numerical comparisons would have to be made between D(T) values to find a global minimum. 3.7 CLINICAL DETERMINATION OF THRESHOLD PARAMETERS In order to use equation (2 5) to optimize the stimulus pulse width to minimize stimulation damage, i t i s necessary to determine numerical values for the threshold parameters a, b and c A method to evaluate these parameters in a c l i n i c a l setting for a biphasic stimulus with given phase delay i s outlined below. The expression for threshold current calculated from the RC-variable threshold assumption i s given by: T —bT _ -bT. I (T ) = — —— + — —. (26) t h (1 - e"1^) (1 - e ) C l i n i c a l l y , one can measure I . and I. directly. It w i l l be shown that 'b' J rh 9 may be readily evaluated from knowledge of the chronaxie, T£. The rheobase c u r r e n t , Ir n» can be e a s i l y obtained by noting the threshold amplitude of long duration (but non-accommodative) pulses. The 1^  parameters must be determined by noting the threshold amplitude of a short duration Dirac delta stimulus of width 6, denoted by 1^(6). Ideally 5 i s 69 Figure 3.4 Pulse Widths Minimizing the Normalized Damage Funct ion Pulse width values minimizing the normalized damage f u n c t i o n , parameterized by phase s e p a r a t i o n , as a f u n c t i o n of the damage parameter. 70 zero, but practically i t is only necessary that i t be a value significantly smaller than the membrane time constant x, (which for the FH model i s about 66 us). With a sufficiently small 6 value, I0 i s defined such that 6 • I . (6) Q (6) I = ^ : (27) 9 T T V ' The i n t e r p r e t a t i o n of I i s that i t represents the charge injected through 9 the membrane as capacitive displacement current in time 6. Mathematically, i t i s related to the instantaneous potential increment required to displace the phase point of the FH system to the 'base' of the threshold separatrix [15]. From the rheobase current one can determine the parameter 'c' using c = I , • R (28) rh The value of 'a' can be calculated from the rheobase and theta currents and the membrane resistance R: a = ( I9 ' Ir h) ' R ( 2 9 ) The decay constant 'b' may be related to chronaxie by (30) b = i - l n c I — I , 9 rh C l i n i c a l l y , the chronaxie T£ i s determined by f i r s t setting the stimulator amplitude to twice rheobase, then by adjusting the stimulus pulse width control u n t i l a setting i s reached for which the chosen amplitude i s the threshold amplitude. Substitution of the stimulator pulse width setting for Tc i n (3) w i l l allow calculation of 'b'. These variable threshold potential 71 parameters may now be used in (25) to evaluate the optimal pulse width once the damage parameter y has been set to the desired value. 3.8 DISCUSSION In the preceeding analysis, i t has been assumed that the threshold response of a myelinated nerve can be modelled by a passive RC network and a modified decaying exponential threshold potential function. Derived from the Frankenhaeuser-Huxley model of myelinated nerve, the potential function gives a threshold current locus which i s i n good agreement with the FH thresholds to a pulse width of about 100 us, and in only f a i r agreement with pulses longer than this to a reasonable limit of about 500 us. Exhibiting the same functional form for monophasic and biphasic pulses, up to a phase separation of about 20 us, the threshold potential function contains three constants, a, b and c, which must be numerically determined for a given phase separation. The dependence of these parameters upon the phase separation of the stimulus i s notable, but the form of this dependence i s not known. Hence any results obtained using this threshold potential function are for some specified value of phase separation, and thus, in the case of minimizing some stimulus damage function, the minimum separation value can only be determined by relative comparisons of the damage function. For equally weighted contributions of power and charge damage compo-nents, damage minima were found to l i e at pulse widths of about 400 us for monophasic stimuli (in f i n i t e phase separation), and 500 us for biphasic pulses with 20 us separation. These pulse separation settings represent upper and lower limits for which the exponential threshold function i s appli-72 cable. Any other given separation between the monophasic and 20 us value w i l l give a minimum between 400 us and 500 us pulse width., for a y of 0.5. This range of values i s just within the 500 us limit beyond which the RC-variable threshold model i s inadequate. Selection of smaller y values which increase the importance of the charge damage component, would give minimum pulse widths at lower values closer to the region of good agreement between the FH and RC-variable threshold models. For a y of 0.1, the monophasic minimum occurs at 157 us, while the biphasic 20 us minimum is near a pulse width of 234 us. A y value of larger than 0.5 would begin to yield long optimal pulse widths past the range of applicability of the RC-variable threshold model, and thus such results should not be considered accurate. The optimum pulse width ranges as functions of y are shown in Fig. 3.4. For biphasic stimuli with phase separations less than 20 us, Butikofer and Lawrence [44] showed Indirectly, using an FH model, that the minimum of charge on a given threshold locus occurred at some non-zero pulse width set-ting. As phase separation decreased below 20 us, a maximum of threshold charge was demonstrated by threshold Dirac delta function stimuli in contrast to the charge minima shown by Dirac pulses with phase separations of 20 us or greater. Such trends lead one to expect minimum damage pulse width l o c i for short separation biphasic pulses to closely resemble those of Fig. 3.4. How-ever, the l o c i would have vertical s h i f t s , so that as y+0 they would inter-sect the ordinate at a non-zero pulse width corresponding to the location of the charge minimum for the stimulus in question. In addition to minimizing stimulus damage due to charge and heat mechanisms, one may also elect to consider optimization of the stimulus with 73 respect to certain design features of the stimulator. One of a number of such features i s the stimulator compliance. A 'cost' function which mini-mizes both damage and the voltage range of current stimulation i s C = YH + (1 " Y)Q + ni(t) (31) Performing the same minimization procedure subject to the same threshold constraints as before yields a 'minimum cost' stimulus waveform which satis-fies the differential equation | i + kL(k2 - \ et / T) I = k3I2 (32) where the are constants, x i s the membrane time constant, and \ is a Lagrange multiplier. This equation i s in the form of a Riccati equation and has a solution t / T -k2t l( t ) = expU-c . e ) - e et / Tr - X T U e V- \ T U -j The prospect of a global solution for this expression allowing determination of the multiplier X i s bleak, but asymptotic solutions for restricted ranges of t may be possible. For pulses of very short duration equation (33) may be simplified, to give: \t I(t) ~ Sg- ( l + (\ - k2)t) (34) 74 3.9 CONCLUSIONS The strength-duration curve derived from the assumption of passive RC neural behaviour in conjunction with a potential threshold function contain-ing three constants, both produces a good match to the numerically determined FH threshold amplitudes and result in a workable analytical definition of neural threshold. The resulting threshold definition i s applicable to both monophasic and biphasic pulses. This approach to the threshold problem allows explicit optimization of pulse widths with respect to a two component damage function which reflects the possibility of power and charge related stimulation damage. The stimulus that minimizes the damage function exhibits the same rising exponential form as the minimum power waveform of Offner. While both waveforms are identical in the limit of in f i n i t e pulse width, they demon-strate greatest dissimilarity for short duration stimuli. In practice such exponential stimuli are not the most convenient to generate; but for pulses of widths of up to one quarter the membrane time constant a constant current stimulus waveform i s a good approximation, and beyond this width, up to a limit of about one membrane time constant, a linearly increasing ramp-like stimulus may be used to approximate the rising exponential. Use of the RC-variable threshold model allows one to calculate an optimal pulse duration from (25) for either monophasic or biphasic waveforms once the threshold parameters a, b, and c have been determined. A proposal has been presented for determining the three parameters i n a c l i n i c a l setting that involves measuring rheobase current, theta current and chronaxie. 75 4. SIMULATION OF HUMAN MEDIAN NERVE THRESHOLD RESPONSE 4.1 INTRODUCTION It i s well known that a stable propagating action potential can be elic i t e d from an axon by electrical stimulation i f the applied stimulus has an amplitude which i s sufficiently large. The minimum stimulus amplitude necessary for generation of an action potential - the threshold amplitude -is a complicated function of many variables. These variables may be cl a s s i -fied according to the nature of their influence upon threshold amplitude. Broadly considered, there are four classes into which the threshold influencing factors may be separated; functional, physical, geometrical, and environmental. The functional class includes the general effect upon thresh-old of stimulus waveshape, and the influence of stimulus parameters such as pulse width, phase separation and frequency. Temperature, local ionic con-centrations and fibre diameter are examples of factors in the physical cate-gory, while the geometrical grouping encompasses such variables as stimulus electrode configuration and the distance of the target axon from the site of stimulation. The environmental class includes variables due to adjacent passive tissues, the presence of additional fibres in the immediate vici n i t y of the target axon, and the presence of membranes acting as barriers to local ionic diffusion. The boundaries of definition of these four groups are necessarily indistinct since parameters of threshold influence belonging to one group may have direct correspondence to influencing factors in another group. 76 If one attempts to simulate the threshold behaviour of an in-vivo axon, variables from a l l four threshold influencing groups must be considered. Several models of the nerve membrane have been proposed which simulate mem-brane dynamics in response to current injection at an Isolated patch of active membrane. Such models, which include the Hodgkin-Huxley (HH) model of unmyelinated nerve [27], the Frankenhaeuser-Huxley (FH) model of myelinated nerve [1], Dodge's myelinated nerve model [28] and others, are best suited to simulations seeking to reveal functional threshold dependencies. Models of this type must be embedded within further models with broader horizons to simulate the influence of factors in the physical, geometrical, and environ-mental threshold groups. In [44], Butlkofer and Lawrence used an FH model to obtain l o c i of threshold amplitude and charge for a biphasic square stimulus, demonstrating the functional variation of threshold with respect to the parameters of the biphasic stimulus. To compare model results with threshold l o c i from human subjects, the same authors devised a psycho-physical experimental protocol to account for environmental and geometrical threshold influences, and derived a uniform scale factor to deal with the effects of temperature. Using simplified FH models to simulate neural activity at eleven sequen-t i a l nodes of Ranvier, McNeal [16] produced a threshold amplitude locus reflecting functional threshold variations for a monophasic constant current stimulus. He also calculated threshold dependence upon fibre diameter. The model proposed by McNeal also allowed determination of some geometrical factors, accounting for both the stimulus electrode geometry, and introducing distance from stimulus site to target axon as a threshold parameter. In 77 [4], Mortimer employed McNeal's model to determine the threshold dependence of myelinated nerve as a function of axon diameter and position with respect to the stimulating electrode. By extending the FH model to encompass more threshold factors, and by appropriate design of an experimental protocol to allow for the effects of geometrical and environmental threshold influences, this paper attempts to demonstrate the a b i l i t y of an enhanced FH model to predict experimentally measured in-vivo thresholds from human median nerve. In the following text, we w i l l f i r s t address some of the general extensions to the FH model which are required, and then discuss specific model enhancements and experimental protocol with respect to percutaneous and transcutaneous stimulation of human median nerve. 4.2 ENHANCEMENT AND EXTENSION OF THE FH MODEL The FH model i s a system of five non-linear ordinary differential equa-tions describing the space-clamped response of a myelinated toad (Xenopus Laevis) nerve to current stimuli injected at a single node of Ranvier. The FH system has recently been shown to be a member of the 'quasi-threshold' class of excitable systems [15]. The model equations are summarized by Frankenhaeuser and Huxley in [1], and by McNeal in [16], and w i l l not be repeated here. Although derived from experimental observation of nerve from the African toad, there i s some evidence that the model membrane emulates the membrane dynamics of human and other nerve fibres over certain stimulus 78 parameter ranges [44], [45]. The enhancements and extensions discussed here are for the purposes of improving the applicability of the model to the human case, and for expanding upon the model range. 4.2.1 Temperature The effects of temperature are imparted to the FH model through a set of temperature correction factors (applied to a l l of the model rate variables and ionic permeabilities) as well as through the temperature dependence of the Goldman factors in expressions for sodium, potassium and unspecific currents. In [40], Frankenhaeuser and Moore published the results of experiments determining the behaviour of sodium and potassium rate variables and permea-bilities as functions of temperature. The temperature dependence was expressed as a set of multiplicative exponential functions of the form • (I) - Q S " 1 0 " 1 0 (1) where T is in degrees Celsius. The temperature functions are dimensionless, and have the effect of changing by factors of Q^Q for every ten degree tem-perature deviation relative to 20°C In [40], Q1Q values for the sodium and potassium systems were presented, with subsequent publication by Franken-haeuser in [53] of Q,0 data for the unspecific current system. The 79 unspecific Q1Q data were not experimentally determined, but were selected by Frankenhaeuser using unstated c r i t e r i a . The temperature dependence exhibited by the three Goldman factors for sodium, potassium and unspecific currents, occurs as a product of a linear temperature component and an exponentially occurring temperature term. With-in the Goldman factors, temperature is expressed in Kelvin degrees and thus the relative change in factor magnitude for a given temperature increment i s much less than the corresponding change undergone by the Celsius based Q1Q factors• The f i r s t extension to the FH system arises because the experimental situation we wish to emulate exists at a temperature significantly higher than the base 20°C of the FH model. Experimental stimulation was at the fourth finger and at the wrist. Both sites were found to have temperatures near 30°C. Such a temperture is a significant extrapolation beyond the maxi-mum 22.5°C of Frankenhaeuser's original data base, but the 33°C value was nevertheless used for a l l predictive model calculations in view of the inavailability of higher temperature data and the great experimental d i f f i -culty in obtaining i t . At 33°C, the normally constant -70 mV (20°C) resting potential w i l l change appreciably due to the resulting change in thermodynamic equilibria regulating the resting potential. Hence, the resting potential must be expressed as a function of temperature in order to calculate the new 33°C value. In order to do t h i s , one must introduce three constants, not origin-al l y given by Frankenhaeuser, which are related to a chloride ion system. These constants, external chloride ion concentration, internal chloride ion 80 concentration and membrane permeability to chloride ions, may be inferred from data quoted by Frankenhaeuser and from well known physiological con-stants. The function using these constants which expresses resting potential as a function of temperature arises from the assumption of a constant field membrane and gives: V ( T ) = R*T l n (Pp(T> + + VT>[K+]Q + PCl< T>[ C 1~li ( 2 ) r F n (P(T> + PN a( T ) ) tN a +] i + PK( T )tK +] i + PC l( T )tC 1" ] o where R is the gas constant, F is Faraday's constant, T is the absolute temperature, the Pj(t) a r e temperature dependent ionic permeabilities (with T in degrees Celsius), the [•] a r e outside membrane ion concentrations, and the [*]^ are the inside ion concentrations. With the exception of the chloride ion permeability and concentrations, a l l other factors are known from measurements made by Frankenhaeuser. If the three chloride constants are not introduced, then i t is not possible to reconcile the resting poten-tial calculated at 20°C using Frankenhaeuser's data, with the value of -70 mV measured at 20°C by Frankenhaeuser. The external chloride ion concentration may be calculated as 110.5 mM from the data given by Frankenhaeuser for his frog Ringer's solution [39]. The internal chloride concentration can be calculated from the chloride Nernst potential at 20°C [62], and is found to be 29.7 mM. The final unknown - chloride ion permeability - is calculated from the measurement by Franken-haeuser of a 20°C resting potential of -70 mV as 1.2(10~7) cm/s. Using 81 equation (2) with T = 33°C, the resting potential i s found to be -75.14 mV. Such a value can significantly change the threshold behaviour of the model from that at 20°C with V = -70 mV. r 4.2.2 Passive Tissues and Adjacent Nodes In order to model r e a l i s t i c c l i n i c a l situations where the stimulating electrodes are far removed from the nodal membrane, one must consider other extensions to the FH model. Such extensions were f i r s t suggested by McNeal [16], when he proposed a model of myelinated nerve which accounted for a multi-node fibre and sur-rounding passive tissues. Using a time invariant RC combination to model the nodal membrane and assuming a resistive axonal core with myelin sheaths of in f i n i t e resistance and zero capacitance, McNeal formed an eleven node fibre and conceptually placed this fibre into an i n f i n i t e , isotropic external medium. He was able to calculate the responses of any of the eleven nodes to current stimuli by simultaneous solution of eleven f i r s t order ordinary d i f -ferential equations representing the evolving transmembrane potentials of the eleven nodes. We implemented a substantially simpler version of McNeal's eleven node fibr e , simulating the presence of only two adjacent passive nodes in an i n f i n i t e , isotropic medium (Fig. 4.1). However, such a simplification enabled us to retain the ful l y active FH system of equations at the central, active, node thus allowing the nodal impedance to vary during the course of a 82 Figure 4.1 Equivalent Circuit of the Three-Node FH Model Equivalent c i r c u i t for the three node fibre in an i n f i n i t e , Isotropic medium. Only the central FH node i s f u l l y active, the adjacent nodes are passive with impedance values taken from the quiescent state of the FH model. 83 current stimulus. The distorting effect of such a fibre upon the applied stimuli i s illustrated for two stimulus waveforms in Fig. 4.2. Illustrated in Fig. 4.2a is a current waveform entering the active FH node resulting from stimulation with a biphasic square waveform. The effects of the changing active node impedance can be clearly seen in Fig. 4.2a. The leading positive phase of the nodal stimulus shows a characteristic RC decay in i t s amplitude. In the case illustrated in Fig. 4.2a the leading stimulus phase i s sufficient to cause an action potential with a concurrent decrease in membrane imped-ance. The effects of this reduced impedance on current entering the active node can be seen during the t r a i l i n g phase of the stimulus in Fig. 2a. The current waveform in Fig. 4.2b is that which enters the active node as a result of application of a ramped biphasic stimulus similar to the stimulus illustrated in Fig. 4.3. This model extension was used to simulate c l i n i c a l situations where the approximation of stimulating electrodes in direct contact with the nodal membrane was not adequate. From the tissue distorted transmembrane current waveforms of Fig. 4.2, i t can be seen that the applied stimulus and the actual stimulating waveform can differ greatly. In [16], McNeal was the f i r s t to demonstrate this fact, although his distorted stimulus waveshapes show none of the v a r i a b i l i t y due to changing membrane impedances because of his constant impedance simplification. 8 8 • •roo 7!u 7M o'.ia. 7!a o.eo T i M (•illi-Mconds) Figure 4.2 Distortion of Injected Current Waveforms The affects of the three node model and the central node's non-constant impedance on the current waveform injected at the central node. a) Current crossing the central node of Ranvier, originating from a biphasic square current stimulus. Note the affects of changing membrane Impedance (time constant) on the two phases of the stimulus. b) Current crossing the central node of Ranvier, originating from a biphasic ramped current stimulus (illustrated in F i g . 2). The actual injected waveform is significantly distorted from i t s form as i t leaves the stimulator. 85 8 o5o oToJ oToe O T I S oTii 7.20 T i n e ( a i m - s e c o n d s ) Figure 4.3 DISA 15E07 Stimulus Waveform A t y p i c a l stimulus waveform o r i g i n a t i n g from two DISA 15E07 constant current s t i m u l a t o r u n i t s with in ter -phase delay provided by a DISA 15G01 delay u n i t . The r i s e time i s approximately 0.8 us per mA/cm 2 . 86 4.3 EXPERIMENTAL DESIGN Due to several essential differences between the neural model and the in-vivo experimental situation - differences which cannot be reconciled with simple model enhancements - judicious definition of the experimental protocol is required. Two of the most significant differences between model and experiment are; the simulation of a space-clamped fibre response compared to an in-vivo saltatory response; an in-vivo stimulation of many nerve fibres and subsequent measurement of a compound action potential (CAP) compared to simulation of a single fibre response. The problem of comparing single fibre model responses to in-vivo bulk action potential recordings can be addressed by determining in-vivo thresh-olds relative to an arbitrarily chosen set point. The CAP of some nerve trunk, measured i n response to a stimulus of amplitude I with some set of pulse parameters p^, can be considered to be a linear superposition of N single fibre responses [63]: y ( l , P _ , t ) = I g (t) (3) i=l where g^(t) i s the contributed response of the it b active f i b r e , and N i s the number of active fibres. Clearly N w i l l depend upon many factors from a l l four threshold groups including the distance from stimulus electrode to fibre bundle, orientation of the bundle, diameter distribution of fibres within the 87 bundle and, of course, the stimulus amplitude and parameter set. If the location and placement of the stimulus electrodes i s not changed, and i f the fibre population remains the same, then any changes in the CAP due stimula-tion with amplitude I^-AI and the same parameter set p^, i s due to loss of the response contribution of the N fibre ( i f AI i s small enough) with: N-1 y ( l -AI,p ,t) = I g (t) < y ( l ,p ,t) (4) i=l Thus the stimulus amplitude Ir i s defined to be the threshold amplitude of the Ntb f i b r e i n the bundle i n response to the stimulus with parameter set Pr. Precisely which f i b r e of the t o t a l bundle population i s so close to threshold that some small amplitude decrement AI results in a subthreshold response, i s determined by the I n i t i a l placement of the stimulating electrodes, the i n i t i a l amplitude Ir, and the pulse parameter set p^. If y(lr> Pr» t ) i s defined to be the reference CAP response, then the threshold of t h — the N f i b r e for a test stimulus with a different parameter set p^ may be obtained by adjusting the test stimulus amplitude It so that: y(It,Pt,t) = y( lr, pr, t ) (5) The r e s u l t i n g threshold amplitude at pt w i l l be defined as the Ntb fibre threshold relative to I at p , and w i l l be given as: r r 88 W V = (6) Inherent in this definition of relative threshold are the assumptions that any two equivalent CAPs correspond to identical firing populations, with identical orders of fibre recruitment. These two factors were addressed by Janko and Trontelj in an excellent study [10]. Using microelectrodes inserted into human ulnar, median and radial nerves, the authors recorded responses to various transcutaneously applied stimuli subsequently showing that when two stimuli were adjusted to the same neural response level, the same active population of fibres resulted with an identical sequence of recruitment. It has been shown by Cooley and Dodge [46] that differences in threshold loci between space-clamped and freely conducting axons can be accounted for, without significant error, by a uniform scale factor. The influence of any uniform scaling factor is eliminated by means of the definition of relative threshold which entails normalization of the threshold amplitude with respect to the set point threshold 1^. The above protocol definition results in the following sequence of steps which must be performed for either percutaneously or transcutaneously applied stimuli: 2. 1. Select a reference parameter set pr» Select a reference stimulus amplitude I . 3. Stimulate with I at p ; record and retain the resulting set point CAP as the reference response. 4. For various test parameter sets, pfc, determine the test amplitude Ifc required to match the recorded reference CAP and calculate I (p ) = I /I . th t t r 89 Three typical near nerve recordings which ill u s t r a t e this procedure are shown in Fig. 4.4, along with a characteristic DISA current stimulus wave-form. The central recording is that obtained from a reference stimulus of 50 us phase duration with 50 us phase separation (50-50-50 setting, p = (50,50,50)). The reference amplitude of Ir(50,50,50) = 25 mA i s chosen using a protocol outlined in the experimental section of this paper. The resulting A-B peak of 105.0 mV i s retained for reference. The upper record-ing represents that due to a test stimulus with pt = (100,20,100). The measured A-p peak which results from application of this stimulus i s 101.3 mV in amplitude, and i s lower than the 105.0 mV reference peak. Hence the 19 mA amplitude of the pulse must be increased to obtain an A-B response equivalent to that of the reference stimulus. Similarly, the lower recording shows the response to a test stimulus of 25 mA at a pulse parameter set of Pt = (20,80,20), but i n this case the test response - and hence the test amplitude - is too large and must be reduced to yield the correct relative threshold amplitude. 4.4 EXPERIMENTS Threshold data was collected over a period of eight months in six separate t r i a l s involving either a percutaneous or transcutaneous stimulation scenario. The same healthy subject was used for a l l of the experimental sessions. 90 Near Nerve Recording Near nerve recordings demonstrating the concept of a reference stimulus (middle recording) and low and high stimulus amplitude s e t t i n g s at pulse parameter sets other than that of the 50-50-50 reference s e t . A t y p i c a l b i p h a s i c ramp stimulus i s shown below the three CAP r e c o r d i n g s . 91 4.4.1 Percutaneous Stimulation 4.4.1.1 Selection of Model With appropriate assumptions concerning stimulus electrodes, fibre geometry and relative magnitude of certain fibre impedances, a single FH node may be used to emulate the response of median nerve fibres to percutaneous stimulation. We f i r s t assume that the fibres of the median nerve are circular in cross section, that they are aligned so that a l l fibre axes are p a r a l l e l , that they are arranged such that a transverse cross-section reveals a circu-lar closest packing geometry (Fig. 4.5), and that a l l fibres are of the same diameter. In addition, we assume that the stimulus current vectors are preferentially oriented transverse to the fibre axes with no longitudinal path component, and that the current vectors are contained within the bounds of a cylindrical element of unit cross-sectional area whose axis i s perpen-dicular to the axis of the fibres being stimulated. With these assumptions, a single FH node may be used to model the percutaneous stimulation scenario. 4.4.1.2 Electrophysiological Techniques Using two DISA model 15E07 constant current stimulators and a DISA model 15G01 delay unit, biphasic stimuli were generated and applied through needle electrodes inserted at the wrist. The stimulus electrodes were positioned as close as possible on either side of the median nerve using anatomical land 92 Figure 4.5 Assumed Geometry f o r Percutaneously A p p l i e d S t i m u l i The assumed geometry, stimulus e lec t rode o r i e n t a t i o n and stimulus current path f o r the case of percutaneously a p p l i e d s t i m u l i . The equivalent c i r c u i t ( a s i n g l e a c t i v e FH node) i s a l so shown. 93 marks. Electrode position was 'fine tuned' by stimulating with a low ampli-tude monophasic pulse while simultaneously adjusting the electrode position until a maximal motor response to the stimulus was attained. Responses to the biphasic stimuli, applied at a rate of 1 Hz, were recorded using a near-nerve technique [64] with electrodes placed near the median nerve at the upper arm. The CAP responses were amplified by a DISA 15C01 amplifier, then passed to a DISA 15G01 signal averager where sixteen responses were averaged to reduce background noise. 4.4.1.3 Stimulus Waveforms It was observed that the stimuli generated by the DISA 15E07 units had non-zero r i s e / f a l l times which were proportional to stimulus amplitude. The measured r i s e / f a l l rate was approximately 0.8 us per milli-amp of current density (Fig. 4.3). The internal circuitry of the 15E07 units was modified to allow investigation of stimulus waveforms with pulse durations below the normal 100 us minimum. With this modification, stimuli with durations ranging from 20 us to 100 us and phase delay times of 20 us to 80 us were investigated. The biphasic pulses were generated by chaining two monophasic pulses with the same duration, but opposite polarity, through the DISA delay unit which was used to control phase separation. 4.4.1.4 Reference Selection The parameter set for the reference stimulus and CAP was chosen to have a 50 us i n i t i a l duration, a 50 us phase separation, and a 50 us duration fi n a l reverse phase. This 50-50-50 parameter set was selected because the 94 values corresponded to approximate mid-ranges of duration and separation for the parameter range under investigation. The reference amplitude was chosen to be the mean of sensory threshold amplitude and A-8 maximum amplitude. The A-8 fibres (largest diameter, fastest conducting myelinated fibres) are known to have the lowest threshold to stimulation [65]. The A-8 maximum stimula-tion amplitude i s the amplitude beyond which the CAP peak corresponding to this fibre population no longer increases. This indicates that there are no further inactive A-8 fibres in the nerve trunk which can be activated by increasing the level of stimulation. Since the A-8 group i s the fastest conducting, i t s peak may be identified as the f i r s t occurring peak in the detected CAP. The reference 50-50-50 stimulus amplitude selected typically lay i n the range of 2.8 mA to 3.4 mA. The reference CAP response to a typi-cal set point stimulus i s shown in Fig. 4.4. 4.4.1.5 Procedure The subject was allowed to rest on a bed in a darkened room, and the stimulus and recording electrodes were then Inserted. Sensory threshold and the A-B maximal stimulus amplitudes were determined and from these values the 50-50-50 reference amplitude was calculated. A reference CAP was then generated, recorded, and stored on one of the four channels of the DISA 15H01 display unit. Relative threshold amplitudes were then calculated for twenty parameter sets corresponding to five pulse widths ranging from 20 us to 100 us at four phase delay settings from 20 us to 80 us. The threshold points were determined by adjusting the test stimulus amplitude until the A-8 peak of the resulting CAP matched the 50-50-50 reference response. As well, 95 stimulus amplitudes required for a just noticeable increase and decrease of the test CAP away from the threshold level were determined- After every five test thresholds the reference CAP was re-evaluated to determine i f there had been a reference change due to, for instance, movement of the stimulus or recording electrodes. If there had been a reference change, the reference CAP was redefined and threshold values relative to the new reference stimulus were determined for the previous five test settings, whereupon the reference was again checked. For the twenty parameter sets involved, this procedure would typically take from two to two and one-half hours. Throughout the c l i n i c a l session, EMG activity was monitored by audio output. In the event that spurious contributions to the CAP were made by EMG signals, that CAP was rejected and another was recorded. If EMG activity persisted and It was not possible to eliminate the activity by encouraging the subject to relax, then the session was terminated. 4.4.1.6 Comparison with the Model The relative threshold amplitudes of model and experimental data are shown in Fig. 4.6. In general, the comparisons between model and experi-mental threshold amplitudes are favourable. The same general trends are exhibited by both data sets, highest threshold at lowest pulse width with shortest phase delay, and lowest threshold at highest pulse width with longest phase delay. The greatest discrepancies occur at low pulse width and .short phase delay, where the model data consistently under-estimates the required threshold current. Allowing that the model thresholds are accurate to within 1% of the true threshold value for each stimulus parameter setting, 96 8 20 us phase separation © 30 us phase separation <• 50 us phase separation BO us phase separation -Model data -Experimental data "VM SJOO 4b.oo ab.oo tb.oo iSo~ Pulse Midth («icro-seconds) oo Ik 00 F i g u r e 4.6 Mode l and C l i n i c a l T h r e s h o l d s f o r P e r c u t a n e o u s S t i m u l i R e l a t i v e a m p l i t u d e s o f p e r c u t a n e o u s l y d e r i v e d model d a t a (dashed l i n e s ) compared to e x p e r i m e n t a l r e l a t i v e t h r e s h o l d s o b t a i n e d f rom pe rcu taneous s t i m u l a t i o n ( s o l i d l i n e s ) . 97 and that the error associated with each of the CAP amplitude readings is 2% -only the model and experimental thresholds at the lowest phase separation value of 20 us f a i l to lay within the bounds of experimental error. The modelling assumption that a l l of the stimulus current is injected at a single node of Ranvier is highly simplified. In practice this will cer-tainly not be the case, some fraction of stimulus current will be directly shunted to the cathode (not passing through the membrane) from the stimulus anode resulting in a more pronounced effect for the shortest duration stimuli. This is demonstrated in Fig. 4.2 where comparisons between pure injection stimulus waveforms and tissue distorted waveforms are shown for biphasic squares and biphasic ramped pulses. The effects of a time varying membrane impedance on trans-membrane cur-rent flow can also be seen in Fig. 4.2. The non-constant membrane impedance (and associated membrane time constant) has the effect of changing the frac-tion of stimultor current diverted through the nodal membrane. This can be seen in Fig. 4.2a where the current entering the active node, resulting from a biphasic square current stimulus, is shown. The leading positive phase of the nodal stimulus shows a characteristic RC amplitude decay. As the node is excited, the resulting decrease in membrane impedance influences the fraction of stimulator current which flows into the node. This changing influence can be detected in the second stimulus phase where an obvious decrease in mem-brane time constant can be observed. The same phenomenon for a ramped biphasic square stimulus is shown in Fig. 4.2b. It can be seen from Fig. 4.2 that in the case of a stimulus with non-zero rise time, the effects of current shunting would be most significant for 98 the briefest pulses. It is possible that the observed differences between the predicted thresholds and the measured thresholds are due to the assumption of complete stimulus injection, the inadequacies of which are most pronounced for the shortest temporal parameters. 4.4.2 Transcutaneous Stimulation 4.4.2.1 Selection of Model Due to the physical arrangement of the transcutaneous stimulating elec-trodes, injected current w i l l pass from cathode to anode in a direction which is preferentially oriented along the longitudinal axis of the median nerve. In doing so, a particular current path w i l l certainly traverse a large number of excitable nodes - this i s in contrast to the percutaneous case where current flows across the longitudinal axis of the median nerve. Hence, for the transcutaneous case, the presence of adjacent nodal membranes was deemed important, and such nearby nodes were included in the model for transcutane-ous stimuli. the nodal geometry, electrode orientation and stimulus current path for this case i s shown in Fig. 4.7. A single active FH node with passive nodes on either side was used to simulate the median nerve response to the transcutaneously applied stimuli. The ele c t r i c a l network representation of this model is shown in Fig. 4.1. Only two nearby nodal membranes were simulated, and they were treated as unexcitable. Both of these simplifications were made to reduce the com-99 Figure 4.7 Assumed Geometry for Transcutaneously Applied Stimuli The assumed geometry, stimulus electrode orientation and stimulus current path for the case of transcutaneously applied stimuli. For c l a r i t y , the stimulus current paths have been omitted. The equivalent c i r c u i t (a single f u l l y active FH node flanked by two passive nodes) Is also shown (this c i r c u i t i s given in more detail in F i g . 4.1). 100 plexity and computation time of the model. The passive nodes were represented by resistive and capacitive values from the resting state of an FH node. 4.4.2.2 Methods and Procedures Only the method of applying the stimulus current was changed for this set of experiments. A l l signal recording and reference selection was identi-cal to the percutaneous stimulation case. Data was collected for the same twenty parameter sets as in the percutaneous experiments, using the same procedure of regularly checking the reference CAP during the threshold measurements. 4.4.2.3 Stimulation Two ring electrodes with surface areas of 4.8 cm2, moistened with a saline solution, were placed around the fourth finger. Stimulus current was applied through these electrodes, and twenty threshold points were measured relative to a 50-50-50 reference setting. As before, EMG activity was monitored with audio output. The reference stimulus amplitude was chosen as the mean of sensory threshold and A-8 maximal amplitudes, with typical set-tings in the range 14.8 mA to 15.4 mA, corresponding to current densities of 3.1 mA/cm2 and 3.2 mA/cm2, respectively. The time taken to obtain the twenty threshold values was generally t«o to two and one-half hours. During the course of the experimental session, a control set of stimulus and recording electrodes was monitored on the subject's other arm for the purposes of detecting any changes in threshold due to drying of the stimulus electrodes. 101 No such changes were observed during the length of time taken by an experi-mental session. 4.4.2.4 Comparison with the Model The relative thresholds of model and experimental data for the transcu-taneous case are shown in Fig. 4.8. As with the percutaneous case, the com-parison of the two threshold data sets i s generally favourable. However, for a l l phase delay settings, the model data at long pulse durations i s consis-tently lower than the measured threshold amplitudes but with only the values at 20 us and 30 us being separated by more than experimental error. As the test pulses become shorter, correspondence of the predicted and measured threshold values is good. This i s in contrast to the general underestimation of the experimental results for short pulses in the percutaneous case. The low predicted threshold amplitudes for long duration pulses (80 us and 100 us duration) are possibly due to the highly simplified network chosen to model the effects of adjacent nodal membranes. As the stimulating pulses are increased in duration, the adjacent nodes begin to respond to the presence of a stimulating current thus influencing the current distribution near the active node. The simple network used accounts only for a single passive node on either side of the active target node. This simplification appears to become less adequate as the applied pulses become longer, and more nodes may have to be simulated in order to improve the accuracy of the model predictions. For pulses of durations up to 1 ms, McNeal [16] determined that ten nodes, five on either side of a target node, were necessary to obtain a highly accurate simulation of fibre a c t i v i t y . It should be noted that 102 8 o •J • c • 20 us phase separation © 30 us phase separation «> 50 us phase separation D BO us phase separation -Model data -Experimental data ^00 80.00 46.00 86.00 86.00 i5o~ Pulse Width (licro-seconds) oo Ik Figure 4.8 Model and C l i n i c a l Thresholds for Transcutaneously Applied Stimuli Relative amplitudes of transcutaneously derived model data (dashed lines) compared to experimental relative thresholds obtained from transcutaneous stimulation. 103 McNeal's simulations took place at a temperature of 20°C, while the data in this simulation correspond to a temperature of 33°C. The pulse widths at which adjacent nodal activity becomes significant w i l l be lower i n a 33°C environment than in a 20°C environment due to an appreciable decrease in model time constants at the higher temperature. 4.5 DISCUSSION With the exception of determining sensory threshold, the subject had no conscious input into the experimental data. It was thus decided that a large number of subjects was not required for this study. Had the determination of threshold amplitude involved a perceptual interpretation of sensation by the subject, then clearly a larger number of subjects would have been required for the study. Over the eight month period during which the data was col-lected, the relative threshold l o c i were found to be quite repeatable. Although the absolute magnitude of the determined thresholds could, and did, vary from session to session, the thresholds relative to the 50-50-50 refer-ence point exhibited l i t t l e v a r i a b i l i t y . During a particular experimental session, variation of the reference amplitude was not large, and more fre-quently than not, there was no variation at a l l . The largest ' d r i f t ' of reference threshold amplitude recorded was 0.4 mA in a transcutaneous ses-sion, and 0.2 mA during a percutaneous session, representing 1.5% and 2% amplitude variations, respectively. These d r i f t s are l i k e l y due to shifts in 104 position of the recording electrodes, or in the percutaneous case, the stimu-lus electrodes as well as the recording electrodes. In the case where assumption of a pure injection waveform was not ade-quate, i t was found that the current waveform crossing the neural membrane was significantly distorted from i t s form as applied through the stimulus electrodes. The distortion proved to be most significant for those stimuli with ramped, non-zero, rise times due to the capacitive properties of the FH nodes. Use of a fu l l y active FH node allowed investigation of the influence of membrane non-linearities upon neural threshold. In previous neural models [16], passive FH nodes have been used to illuminate the effects of stimulus electrode and fibre geometry upon neural threshold. However, passive FH nodes represented by time-invariant RC networks have the disadvantage that their threshold i s independent of the tr a i l i n g phase of any biphasic stimu-lus. This limits the scope of any waveform investigations which may be undertaken with such a model. Earlier applications of the FH model for threshold response prediction in humans [44], have used less fundamental means to account for the influence of temperature upon threshold. While global temperature scaling techniques of the style employed in [44] are adequate, they prevent detailed examination of the effects of temperature on the dynamics of the neural membrane. Only a subset of the f u l l geometric considerations addressed in [16] has been implemented with this FH modelling system, in favour of greater atten-tion to the details of membrane non-linearity and time variance. Although emphasis has been placed upon the dynamics of the FH membrane, geometric 105 factors are s t i l l important as can be seen when comparing the c l i n i c a l results in Fig. 4.6 and Fig. 4.8. These figures represent two fundamentally different c l i n i c a l stimulation geometries, and the c l i n i c a l l y recorded thresholds reveal the underlying geometric influences in both the rate of change of reltive threshold amplitude with respect to pulse width (a lesser negative value in the percutaneous case than in the transcutaneous case), and the rate of change of relative threshold amplitude with respect to phase separation (greater in the transcutaneous case than in the percutaneous case) . A comparison between the two experimental relative threshold l o c i for phase separations of 20 us and 80 us i s shown in Fig. 4.9. The changes made to the FH model greatly enhanced Its a b i l i t y to simu-late the experimental situations under investigation. Although the tempera-ture Q10's determined by Frankenhaeuser are, s t r i c t l y speaking, valid only to 22.5°C, the extension of this data base into temperature ranges more relevant to human nerves was made because the Q0^ values were shown to be constant over a range of two decades in Frankenhaeuser' s data. The simulation of two adjacent passive nodes resulted in a significant increase in model-relative threshold at low stimulus durations, yielding threshold values much closer to those of the transcutaneous scenario. 106 F i g u r e 4.9 C l i n i c a l Thresholds f o r Percutaneously and Transcutaneously A p p l i e d S t i m u l i C o m p a r i s o n of r e l a t i v e t h r e s h o l d a m p l i t u d e s f o r c l i n i c a l l y a p p l i e d transcutaneous and percutaneous ramped constant current s t i m u l i . The threshold l o c i are phase separations of 20 ps and 80 u s . 107 4.6 CONCLUSIONS In order to use the Frankenhaeuser-Huxley model of myelinated nerve to simulate in-vivo human median nerve thresholds, discrepancies between model and experiment relating to four groupings of threshold influencing factors must be addressed. The differences may be accounted for by extending the temperature range of the FH model to encompass typical temperatures of human nerve. The temperature extension is achieved through introduction of a constant chloride system which allows the resting potential to be expressed as a function of temperature. This resulted in a value of Vr = -75.14 at the c l i n i c a l temperature of 33°C. Such a potential change, although not large in magnitude, w i l l strongly alter the threshold characteristics of the FH sys-tem. The mechanism of these threshold changes can be understood by examining the effects of a changed V"r on the FH phase plane [15]. For transcutaneous stimulation scenarios, where direct current injection at a single active node i s not a good model approximation, the model must be extended to include the effects of adjacent nodes and passive tissues. In addition, definition of an experimental protocol which yields threshold data relative to an arbitrarily chosen reference stimulus, eliminates threshold influencing factors such as stimulus and recording electrode placement, and the uniform saltatory conduction scale factor. With the above model enhancements and extensions, transcutaneous and percutaneous threshold data collected from human subjects under the 'set point / reference' protocol, closely matches the thresholds predicted by the FH model for biphasic ramped current pulses. The discrepancies between model 108 and experimental results may be due to stimulus current shunting and the influence of nearby nodal membranes for the percutaneous and transcutaneous cases, respectively. Given the nature of the stimulus waveform distortions caused by adjacent nodes and passive surrounding tissues, stimuli with as l i t t l e ramping as possible are preferred. Square stimuli undergo l i t t l e distortion compared to ramped stimuli, and hence the square pulses are more e f f i c i e n t . Consequently a factor of significance in stimulator design is the rise time parameter of the power transistors in the output stage of the stimulator - this rise time should be minimized. 109 5. SIMULATION OF NODAL POTASSIUM ACCUMULATION 5.1 INTRODUCTION It i s well known that the outward potassium current which occurs during the repolarization phase of an action potential can give rise to an accumula-tion of potassium ions in the immediate v i c i n i t y of the active membrane [66], [67]. Although the phenomenon was i n i t i a l l y noticed in the squid giant axon, the same effect has also been observed in myelinated nerve from frog (Rana esculenta) [68], [69], [70], and in the snail (Helix pomatia) [71]. Differ-ences in the mechanism and parameters of K+ ion accumulation and diffusion have been proposed in [69] as one of the reasons for the observable di f f e r -ences between motor and sensory nerve action potentials. Two accumulation models were proposed by Frankenhaeuser and Hodgkin [66], in the original work on K1" accumulation. Both models are s t i l l accepted, with minor alterations, i n more recent literature [68]. Using the nomenclature of Moran et. a l . [68], the two mechanisms are identified as the 'three compartment model' (TCM) (or multi compartment model), and as the 'diffusion in an unstirred layer' (DUSL) model. The work presented here w i l l be concerned with how these accumulation models simulate K+ ion accumulation i n such a way as to effect the threshold c h a r a c t e r i s t i c s of the Frankenhaeuser-Huxley axon. 110 5.2 MODELS OF K+ ACCUMULATION 5.2.1 The Three Compartment Model Originally proposed by Frankenhaeuser and Hodgkin for the squid giant axon in [66] as 'Hypothesis 1: f i n i t e space and very thin barrier to diffusion' , i t has become known in the most recent literature as the three compartment model of diffusion (TCM). The TCM, which consists of the fi b r e , the external space and the bulk solution, explains K"1" accumulation in terms of diffusion obstruction caused by other anatomical structures (membranes). In this model, the outward flowing K+ ions responsible for the potassium current are blocked by a membrane, some distance 8 from the axon, which has a permeability, PK, to potassium such that ions escape from the external compartment at a rate slower than their rate of entry into the compartment. The result i s that the active membrane 'sees' an inflated external potassium ion concentration, with a subsequently increased K+ Nernst potential. Analysis of their squid giant axon data led Frankenhaeuser and Hodgkin to the conclusion that the TCM was the best description of their experimental situation. They calculated the space thickness, 9, of the external compartment to be 300 Angstroms, with an apparent K+ permeability of 6*10"5 cm/s. In [72], Adelman and Palti reported a space thickness of about 360 Angstroms with a permeability of 3.2*!©"1* cm/s. There is anatomical evidence to support the assumption that the K+ ion accumulation is due to nearby membranes as the 9 values correspond f a i r l y well to known Schwann space dimensions in the squid. I l l 5.2.2 Diffusion in an Unstirred Layer Model This model was put forth by Frankenhaeuser and Hodgkin [66] as 'Hypothesis 2: f i n i t e diffusion barrier and no space', along with their TCM hypothesis. In the DUSL model, an elevation in external potassium ion con-centration is observed because of the relatively slow mixing of the bulk solution with the contents of an aqueous layer adjacent to the active mem-brane patch from which the K+ ions emerge. This model, thermodynamic rather than anatomical in origin, requires two parameters to describe the accumula-tion phenomenon - the unstirred layer depth, and the effective unstirred layer K"1" diffusion constant. 5.3 CHOICE OF A K+ ACCUMULATION MODEL FOR THE FH SYSTEM Based upon the recent findings of Moran et. a l . [68], the DUSL model of accumulation was chosen to simulate JC*" accumulation in the FH model. In their experiments, Moran et. a l . found that data for myelinated frog nerve (Rana esculenta) was well fitted by a TCM - this was determined by plotting' TCM theoretical potassium Nernst potentials, for various voltage clamp depolarizations, against the duration of depolarization, and comparing these values to the experimentally measured potentials. However, the anatomical basis for this model i s not obvious since, for a frog myelinated axon, a continuous structure which surrounds the axon has not been described - as reported by Moran et. a l . Since the TCM had no apparent structural basis, Moran et. a l . postulated that the IT*" accumulation was due '...solely to the discontinuity of ion-transport number in the path of electrical current flow, accompanied by slow 112 mixing at the (nodal) surface' . They tested this hypothesis by comparing potassium Nernst potential shifts as functions of depolarization duration to DUSL calculations, and concluded that for 17 fibres the averaged DUSL parameters which best corresponded to the experimental data were a DUSL thickness of 1.4 um, and a diffusion constant of 1.8*10~6 cm2/s. The DUSL data so calculated was observed to f i t the transient portion of the K+ Nernst potential shift (for up to 10 ms at 15°C), but to exhibit less agreement with the steady state experimental value. This i s in contrast to the observation of Frankenhaeuser and Hodgkin that the DUSL model did not accurately describe K+ accumulation around the squid giant axon for short time periods, but was more accurate in the limit of long depolarizations. The DUSL model was selected for this simulation of accumulation because of the observed agreement with experiment at short (< 10 ms) depolarization times, and because of lack of anatomical evidence for a TCM. 5.4 THE AFFECTS OF K+ DUSL ACCUMULATION ON THE FH SYSTEM The extent of potassium ion accumulation during the 'standard' FH response to a 120 us monophasic current stimulus of 1 mA/cm2 is summarized in Figs. 5.1 and 5.2. In Fig. 5.1, the time course of the outward potassium current i s shown, while Fig. 5.2a details the resulting increase in external potassium ion concentration. The corresponding changes in resting potential are illustrated in Fig. 5.2b. For the calculation of Fig. 5.2, the parameter values of Moran et. a l . are used, as well as parameter values for the DUSL system in which no IC*" ions are allowed to escape from the unstirred layer. Assuming the standard response to be quite typical, the latter curves show 113 the theoretical maximum potassium accumulation, and corresponding maximum effect, for a single action potential. The external potassium ion concentration appears directly i n only one of the standard FH equations - that describing potassium current: (V+Vr) - [K+]o*e(V+Vr)F/RT h = PK( V ) F ~RT~ 1 - e(V-W JF/RT However, as previously mentioned, the constant -70 mV resting potential may be replaced by an expression which allows calculation of at some tempera-ture T based upon permeabilities and ionic concentrations of sodium, potas-sium and chloride ions: V (T,[K+] ) = — *l n ^ P ^ ^ N a ^ ^ t ^ l o ^ ^ ^ o ^ q C ^ t0 1^ ! , r ' ° F (Pp(T)+PN a(T))[Na+]1+PK(T)[K+]1+Pc l(T)[Cl-]o This expression of Vr as a function of external K"1" concentration i s essential i f the f u l l influence of K"1" DUSL accumulation is to be taken into account. The global effects of increased external K+ concentration on the FH system are best observed in the reduced FH phase plane; this i s shown in Fig. 5.3, where the effect of a 20% increase in external K**" concentration upon a mem-brane in i t s resting state is demonstrated. As far as neural threshold is concerned, F i g . 5.3 shows that an increase in w i l l have the effect of raising threshold by virtue of the translocation of the unstable saddle point 114 8 Figure 5.1 Potassium Current Flow in Response to Standard Stimulus Membrane potassium current flow in response to the 'standard* FH stimulus, 1 mA/cm2 120 us duration monophasic constant current sitmulus. 115 2.40 .00 oT*©Teoi;» Tiie (aiHi-seconds) Figure 5.2 FH Variables Affected by Potassium Ion Accumulation a) External potassium ion concentration during FH standard response using DUSL data of Moran et. a l . , and assuming no diffusion from the unstirred layer. b) Resting potential changes during FH standard response using DUSL data of Moran et. a l . , and assuming no diffusion from the unstirred layer. 116 0.08 (volts) Figure 5.3 Phase Plane Changes Due to Potassium Ion Accumulation The quiescent FH Vm phase plane (thick lines) showing phase plane changes brought about by a 20% Increase in external potassium ion accumulation (thin l i n e s ) . 117 (and i t s associated separatrix) towards larger values for V and m. This translocation is shown in greater detail in Fig. 5.4. Further, the unstable excited node of the system is shifted to a smaller V value, implying a reduced maximum amplitude of the action potential. Such amplitude reductions have been observed in action potentials of the squid giant axon [66], [67]. It i s likely that other FH variables, such as some of the sodium i n a c t i v a t i o n rate constants, are influenced by changing [ K+]Q values. However, these dependencies must be illuminated through further voltage clamp studies as they cannot be deduced from existing FH data. Such additional voltage clamp studies have been done for the Hodgkin-Huxley system, and dependencies of the sodium inactivation variable (h) rate parameters, time constant and steady state values have been revealed [68], [73]. Similar alterations to the FH system without benefit of further voltage clamp data could only be made in a speculative manner, and hence the present modifica-tions of the FH system to account for K+ accumulation should be considered as only f i r s t order in extent. 5.5 THE INFLUENCE OF K+ ACCUMULATION ON MULTIPLE PULSE THRESHOLDS The phenomenon of K*- accumulation can only occur when the outward mem-brane potassium current becomes significant. Hence, K+ accumulation w i l l never influence the threshold characteristics of an excitable node being stimulated from i t s quiescent state since current becomes significant only during the latter (repolarization) phase of the action potential. Potassium accumulation w i l l thus only effect neural threshold in a multiple stimulus -multiple response situation, where the accumulted outside a node as a 118 Figure 5.4 Phase Plane Detail Detail of phase plane of Fig. 5.3 showing unstable saddle point 119 result of a primary response may alter the threshold characteristics of the node with respect to t r a i l i n g stimuli. The simplest means by which to test the degree of potassium accumulation influence is to f i r s t condition the FH system by stimulation with a supra-threshold pulse, then to calculte the minimum (threshold) delay time required in order that a t r a i l i n g stimulus, identical to the conditioning pulse, successfully activate the model node a second time. Since accumultion of extra-nodal potassium has the effect of raising the required threshold ampli-tude at any given stimulus parameter set, one would expect to observe increasing threshold delay times as the amplitude of the tr a i l i n g pulse approached the absolute threshold amplitude for the pulse parameter set under consideration. The threshold delay times calculated with simulation of external potassium accumulation would be expected to be larger than the delay times determined during the relative refractory period of a nodal membrane without simulation of K"1" DUSL accumulation. This i s shown in Fig. 5.5, where the threshold delay time difference between FH systems with and without DUSL simulation i s plotted as a function of the percentage threshold amplitude of the conditioning and tr a i l i n g pulse for 100 us monophasic square constant current stimuli. As anticipated, the required minimum delay times decrease as stimulus amplitude increases above threshold. A vertical asymptote appears at a threshold amplitude fraction of 1.0 since any amplitude fraction less than this value is subthreshold and hence can never activate the excit-able node. For pulse amplitudes closer to threshold, the delay difference i s quite pronounced - on the order of 1 ms. As pulse amplitudes become larger, the small threshold Increase due to IC*" accumulation becomes less significant, 120 Figure 5.5 Threshold Delay Times Threshold delay times of tr a i l i n g 100 us duration monophasic constant current stimulus (threshold amplitude 0.8672 mA/cm2) after an identical conditioning stimulus; as a function of threshold amplitude fraction. 121 and the delay difference decreases. For an FH system with no DUSL simula-tion, a t r a i l i n g pulse with an amplitude within 1% of absolute threshold has a threshold delay of 15.51 ms. With DUSL simulation, the delay threshold was 16.406 ms. The largest amplitude pulse tested, about 1.7 times threshold, had a threshold delay time of only 1. 798 ms, with DUSL simulation yielding a 1.803 ms delay. 5.6 DISCUSSION The simulation chosen to demonstrate the effect of potassium accumula-tion has no experimental twin for which data has been published. The usual method of experimental evaluation of potassium accumulation is to stimulate over extremely long time intervals, and to observe the gradual changes brought about by potassium accumulation. The stimulation times for such experiments are on the order of hundreds of milli-seconds, and hence model emulation of these results i s not practical as calculation of a single 200 ms model response would consume approximately 10 hours of computer time on a PDP 11/23. In [66], Frankenhaeuser and Hodgkin performed a set of experiments which are quite amenable to simulation. They measured the change in t r a i l i n g action potential positive phase amplitude at varying time intervals following a conditioning stimulus. Unfortunately, similar experiments have not been reported using myelinated fibres. The modelling of the potassium phenomenon reported here should be con-sidered only as a preliminary attempt to account for the influence of ionic accumultion upon the FH system. In as much as this phenomenon produces a 122 long term multiple-pulse/multiple-response influence on the FH system, i t represents a departure from the shorter term phenomena with which this work is principally concerned. Much more model modification in conjunction with further experimentation is necessary to discover what additional model parameters have [ K+]Q dependencies. Until such work has been accomplished, these results, while serving to indicate the trends of influence of the potassium assumulation extension to the FH model, do not represent absolute data intended for direct comparison with parallel experimental work. 123 6. SUMMARY This thesis has studied the threshold behaviour of myelinated nerve by using a computational model, the Frankenhaeuser-Huxley (FH) model. The ques-tion of applicability of the model results to threshold properties of human peripheral nerve was also addressed. Chapter 2 began by examining the mathematical properties of the FH system of differential equations in terms of st a b i l i t y and ex c i t a b i l i t y , with particular reference to the nature of threshold in the system. This examina-tion of the FH system used techniques of phase-space analysis, and suggested the form of a computational algorithm for use in threshold calculations on a dig i t a l computer. With knowledge of the FH system's threshold properties, i t was possible to suggest methods to improve the modelling of nerve excitation in a way which accounted for intrinsic membrane dynamics, but yet which was not excessively complicated. The performance of the model was then evaluated against c l i n i c a l l y derived thresholds from human median nerve. To effect the comparison with c l i n i c a l data, i t was necessary to enhance the computational capabilities of the FH system, and to introduce a protocol for collecting in-vivo human threshold data. The comparison of model to c l i n i c a l data was found to be favourable, with the model predictions lying within experimental error of the c l i n i c a l l y measured thresholds. Finally, a preliminary mechanism was provided for enhancing the computational capabilities of the FH system to include the affects of extra-nodal potassium ion accumulation. Phase-space analysis of the FH system revealed i t to be a member of the quasi-threshold class of excitable systems for physically meaningful tempera-124 tures. Using the sharpness definition of Fitzhugh to quantify the degree of excitability sharpness, i t was possible to show that the FH system is two orders of magnitude 'sharper' than another familiar nerve model - the Hodgkin-Huxley model of unmyelinated nerve. The sharpness of the FH system was shown to be such that no gradation of stimulus-response curves was observed up to temperatures of 40°C Analysis of a reduced FH phase space, the plane of the variables V (transmembrane potential) and m (sodium activa-tion variable), revealed that the excitability of the FH system was dependent upon the location of a saddle point and i t s associated separatrix in the Vm phase plane. A least-squares surface was fitted to the separatrix as a func-tion of the i n i t i a l l y most active FH variables V, m and h (the sodium inactivation variable). This surface was used to automate a separatrix threshold definition, which was found to be unambiguous in cases where previ-ous definitions, the level and inflection definitions, f a i l e d . Computer implementation of the separatrix definition was found to significantly reduce computation time when iterating to a value of threshold stimulus amplitude. An improvement on the modelling of nerve excitation was suggested by the nature of the excitability of the FH system. While previous excitation models assumed a constant potential threshold, the FH system indicated that potential thresholds were not constant, but rather were strong functions of stimulus parameter sets. It was shown that the FH variable thresholds could be represented by exponential curves, with three parameters, for biphasic and monophasic square constant current stimuli. The general expression for the variable potential threshold was used in conjunction with a simplified passive RC model of nerve to produce analytic expressions for threshold 125 current amplitudes of square stimuli as functions of stimulus duration and inter-phase delay. It was then shown how such expressions could be used to minimize various 'cost' or 'damage' functions, which yielded expressions for minimum charge and amplitude-squared (proportional to ohmic heating) damage components of varying proportion. The optimizing procedure resulted in an analytic expression for an exponentially rising stimulus wave shape with a pulse duration setting parameterized by phase separation for various 'weight-ings' of the charge and heating damage components. This expression was solved to produce a chart of optimal pulse widths. How the derived optimiz-ing procedure could be applied to optimize any cost function was then shown, and as an example, the optimization of a cost function containing components for charge and heating damage as well as a stimulator compliance components was demonstrated. The c l i n i c a l relevance of this mathematical procedure was indicated by i l l u s t r a t i n g how to c l i n i c a l l y determine the three threshold parameters used in the optimizing process. The three parameters were shown to be related to the familiar rheobase current, the chronaxie, and to a 'theta current' which was introduced and defined in Chapter 5. The theta current was shown to bear a relation to the charge injected through the nerve membrane, as capacitive displacement current, in an element of time which was short in comparison with the membrane time constant. C l i n i c a l l y , i t was demonstrated that the FH model could predict the threshold functions of in-vivo human median nerve, within experimental error, in response to ramped current stimuli generated by a DISA Neuro-Myograph. The stimuli were applied both transcutaneously and percutaneously, and the resulting compound action potentials were measured using a near-nerve record-126 ing technique. In order to compare the model to c l i n i c a l threshold data, i t was shown that threshold-influencing factors from four groups (functional, physical, geometrical and environmental) had to be taken into account. This was done by enhancing and extending the computational capabilities of the FH model, and by defining a set-point/reference (SPR) protocol under which the c l i n i c a l data was collected. The SPR protocol assumed that active axons generate potential fields which are linearly superimposed to form the c l i n i -cally measured compound action potential. Further assumptions regarding recruitment order are required to yield the f i n a l protocol which enables c l i n i c a l thresholds to be compared to model thresholds relative to an arbi-t r a r i l y chosen reference stimulus. Extensions to the model allowed the resting potential to be expressed as a function of temperature and ionic concentrations, and enabled the effects of adjacent passive nodes to be taken into account. The c l i n i c a l and model threshold data were compared for ramped current stimuli with durations from 20 us to 100 us and inter-phase delays from 20 us to 80 us. For percutaneous stimuli, the greatest discrepancies between c l i n i c a l data and model predictions occurred at low pulse widths with low phase separation. These discrepancies are thought to arise from an assumption of complete stimulus injection into the active node. In the transcutaneous case, stimuli with large durations and large phase separations were found to be the only settings where model and c l i n i c a l data were not compatible within experimental error. It Is possible that neglecting the presence of a l l but two adjacent excitable nodes results in the observed long duration discrepancies. 127 Chapter 5 showed how to modify the FH equations to account for the accumulation of potassium ions in the extra-nodal space. The modification involved expression of the transmembrane potential as a function of the internal and external concentrations of three ionic species, K"1", Na+ and Cl~. A thermodynamic model for IC*" accumulation, the 'diffusion in an unstirred layer' (DUSL) model, was incorporated into the FH system and was used to model the affects of accumulation. The affects of K+ accumulation for a single response were shown to be in terms of both increased external IC*" con-centration and increased resting potential. The implications of these affects with regard to model excitability and stability were shown with phase-space analysis. Increasing threshold was shown to result from K+ accumulation, and was demonstrated by measuring the threshold of a monophasic square stimulus which was applied after a supra-threshold conditioning stimulus• In summary, the main contributions of this thesis are: 1. Derivation of the non-constant threshold potential equation, and use of this equation to optimize neural stimuli with respect to a 'cost' func-tion (Chapter 3). 2. Definition of the set-point/reference experimental protocol (Chapter 4). 3. Demonstration of the correspondence of FH model threshold predictions to human nerve thresholds (Chapter 4). 4. Demonstration of the quasi-threshold nature of the FH system (Chapter 2). 5. Derivation of extensions to the FH system involving the potassium system, and resting potential (Chapters 4 and 5). 128 REFERENCES 1. B. Frankenhaeuser and A.F. Huxley, "The action potential in the myeli-nated nerve fibre of Xenopus laevis as computer on the basis of voltage clamp data", J . Physiol. (London), Vol. 171, pp. 302-315, 1962. 2. M.L. Barr and J.A. Kiernan, The Human Nervous System, Fourth Edition, Harper and Row, Philadelphia, 1983. 3. K.B. Roberts, P.D. Lawrence and A. Eisen, "Dispersion of the somato-sensory evoked potential (SEP) in multiple sclerosis", IEEE Trans. on Biomed. Engg., Vol. BME-30, No. 6, pp. 360-364, 1983. 4. J.T. Mortimer, "Motor prostheses", i n : "Handbook of physiology - the nervous system I I " , Vol. 2, Part 1, Vernon Brooks Ed., pp. 155-187, 1981. 5. J. Holle, H. Thoma, M. Frey, H. Bruber, H. Kern and G. Schwanda, "Epi-neural electrode implantation for el e c t r i c a l l y induced mobilization of paraplegics", Proc. 1st Vienna I n t l . Conf. on FES, paper 5.3, October 1983. 6. R. Butikofer and P.D. Lawrence, "Electrocutaneous nerve stimulation II: Model and experiment", IEEE Trans, on Biomed. Engg., Vol. BME-26, No. 2, pp. 69-75, 1979. 7. G. Vossius, "Application of functional muscle stimulation to the paralysed handicapped", Proc. 1st I n t l . Vienna Conf. on FES, Section V, 1983. 8. J. McHardy, D. Geller and S.B. Brummer, "An approach to corrosion control during electrical stimulation", Ann. Biomed. Engg., Vol. 5, pp. 144-149, 1977. 9. R.H. Pudenz, L.A. Bullara, S. Jacques and F.T. Hambrecht, "Electrical stimulation of the brain III: The neural damage model", Surg. Neurol., Vol. 4, pp. 265-270, 1975. 10. M. Janko and J.V. Trontelj, "Transcutaneous electrical nerve stimula-tion: A microneurographic and perceptual study", Pain, Vol. 9, pp. 219-230, 1980. 11. C. Burton and D.D. Maurer, "Pain suppression by transcutaneous elec-tronic stimulation", IEEE Trans, on Biomed. Engg., Vol. BME-21, No. 2, pp. 81-88 , 1974. 12. J.L. Mason and N.A.M. MacKay, "Pain sensations associated with electro-cutaneous stimulation", IEEE Trans, on Biomed. Engg., Vol. BME-23, No. 5, pp. 405-409, 1976. 129 13. R.H. Pudenz, L.A. Bullara, D. Dru and A. Talalla, "Electrical stimula-tion of the brain II: Effects of the blood-brain barrier", Surg. Neurol., Vol. 4, pp. 265-270, 1975. 14. D. Dean and P.D. Lawrence, "Simulation of human median nerve threshold response", submitted for publication, IEEE Trans, on Biomed. Engg. 15. D. Dean and P.D. Lawrence, "Application of phase analysis of the Frankenhaeuser-Huxley equations to determine threshold stimulus ampli-tudes", IEEE Trans, on Biomed. Engg., Vol. BME-30, pp. 810-819, Dec. 1983. 16. D.R. McNeal, "Analysis of a model for excitation of myelinated nerve", IEEE Trans, on Biomed. Engg., Vol. BME-23, No. 4, pp. 329-337, 1976. 17. P.H. Gorman and J.T. Mortimer, "The effect of stimulus parameters on the recruitment characteristics of direct nerve stimulation", IEEE Trans. Biomed. Engg., Vol. BME-30, No. 7, pp. 407-414, 1983. 18. D. Dean and P.D. Lawrence, "Optimization of neural stimuli based upon a variable threshold potential", in press, IEEE Trans, on Biomed. Engg. 19. L. Vodovnik, T. Bajd, A. Kralj, F. Gracanin and P. Strojnik, "Functional electrical stimulation for control of locomotor systems", CRC Critical Reviews in Bioengineering, pp. 63-131, 1981. 20. U. Stanic, R. Acimovic, N. Gros, A. Trnkoczy, T. Bajd and M. Klajajic, "Multichannel electrical stimulation for correction of hemiplegic gait", Scand. J. Rehab. Med., Vol. 10, pp. 75-92, 1978. 21. A. Kralj, T. Bajd, R. Turk, J. Krajnik and H. Benko, "Gail restoration in paraplegic patients; a feasibility demonstration using multichannel surface FES", J. Rehab. R.&D., Vol. 20, pp. 3-20, 1983. 22. H. Stohr, M. Frey, J. Holle, H. Kern, G. Schwenda and H. Thoma, "Func-tional electrostimulation makes paraplegic patients walk again", Proc. 1st Vienna Intl. Conf. on FES, paper 5.4, Oct. 1983. 23. F. Gracanin, "Functional electrical stimulation in control of motor output and movements", Contemp. Clin. Neurophysiol, (EEG Suppl. 34), W.A. Cobb, H.V. Van Duijn Eds., Elselvier, Amsterdam, pp. 355-368, 1978. 24. S. Rebersek, L. Vodovnik, N. Gros, A. St ephanovska, R. Acimovic, "Electrotherapy of spastic muscles in hemiplegia", Advances Ext. Control Human Extremities, ETAN, Belgrade, pp. 217-228, 1981. 25. B.R. Bowman and T. Bajd, "Influence of electrical stimulation on skeletal muscle spasticity", Proc. Adv. Ext. Control Human Extremities, Dubrovnik, pp. 567-576, 1981. 130 26. K. S. Cole, Membranes, Ions and Impulses, Univ. of California Press, Berkley, 1968. 27. A.L. Hodgkin and A.F. Huxley, "A quantitative description of membrane current and its application to conduction and excitation in nerve", J. Physiol. (London), Vol. 117, pp. 500-512, 1952. 28. F.A. Dodge, A Study of Ionic Permeability Changes Underlying Excitation  in Myelinated Nerve Fibres of the Frog, Thesis, The Rockefeller Insti-tute, New York, Univ. Microfilms, Ann Arbor Michigan, (No. 64-7333), 1963. 29. D. Noble, "A modification of the Hodgkin-Huxley equations applicable to Purkinje fibre action and pace-maker potentials", J. Physiol. (London), Vol. 160, pp. 317-352, 1962. 30. R. Fitzhugh, "Dimensional analysis of nerve models", J. Theor. Biol., Vol. 40, pp. 517-541, 1973. 31. L. Goldman and J.S. Albus, "Computation of impulse conduction in myelinated fibres; theoretical basis of the velocity-diameter relation", Biophysical J., Vol. 8, pp. 596-607, 1968. 32. F. Buchthal and A. Rosenfalck, "Evoked action potentials and conduction velocity in human sensory nerves", Brain Res., Vol. 3, pp. 1-122, 1966. 33. F.A. Dodge and B. Frankenhaeuser, "Membrane currents in isolated frog nerve fibre under voltage clamp conditions," J. Physiol., vol. 143, pp. 76-90, 1958. 34. F.A. Dodge and B. Frankenhaeuser, "Sodium currents in the myelinated nerve fibre of Xenopus Laevis investigated with the voltage clamp tech-nique," J. Physiol., vol. 148, pp. ,188-200, 1959. 35. B. Frankenhaeuser, "Quantitative description of sodium currents in mye-linated nerve fibres of Xenopus Laevis," J. Physiol., vol. 151, pp. 491-501, 1960. 36. B. Frankenhaeuser, "Delayed currents in myelinated nerve fibres of Xenopus Laevis investigated with voltage clamp technique," J. Physiol., vol. 160, pp. 40-45, 1962. 37. B. Frankenhaeuser, "Instantaneous potassium currents in myelinated nerve fibres of Xenopus Laevis," J. Physiol., vol. 160, pp. 46-53, 1962. 38. B. Frankenhaeuser, "Potassium permeability in myelinated nerve fibres of Xenopus Laevis," J. Physiol., vol. 160, pp. 54-61, 1962. 39. B. Frankenhaeuser, "A quantitative description of potassium currents in myelinated nerve fibres of Xenopus Laevis," J. Physiol., vol. 169, pp. 455-451, 1963. 131 40. B. Frankenhaeuser and L.E. Moore, "The effect of temperature on the sodium and potassium permeability changes in myelinated nerve fibres of Xenopus Laevis," J. Physiol., vol. 169, pp. 431-437, 1963. 41. B. Frankenhaeuser and L.E. Moore, "The specificity of the in i t i a l current in myelinated nerve fibres of Xenopus Laevis. Voltage clamp experiments," J. Physiol., vol. 169, pp. 438-444, 1963. 42. B. Frankenhaeuser, "Inactivation of the sodium carrying mechanism in myelinated nerve fibres of Xenopus Laevis," J. Physiol., vol. 169, pp. 445-451, 1963. 43. D.E. Goldman, "Potential, impedance, and rectification in membranes," J. Gen. Physiol., vol. 26, pp. 37-59, 1943. 44. R. Butikofer and P.D. Lawrence, "Electrocutaneous nerve stimulation-1: . Model and experiment," IEEE Trans, on Biomed. Eng., vol. BME-25, no. 6, pp. 526-531, 1978. 45. B.M. Clopton, F.A. Spelman, I. Glass, B.E. Pfingst, J.M. Miller, P.D. Lawrence, and D. Dean, "Neural encoding of electrical signals." Pre-sented at the International Cochlear Implant Conference, New York Academy of Sciences, New York City, April 14-16, 1982. 46. J.W. Cooley and F.A. Dodge Jr., "Digital computer solutions for edxcita-tion and propagation of the nerve impulse," Biophysical J. , vol. 216, pp. 932-938, 1969. 47. R. Butikofer, "Electrocutaneous stimulation via bipolar current pulsesL models and experiments," M.A.Sc. Thesis, Appendix A, pp. A - l l . 48. R. FitzHugh, "Mathematical models of threshold phenomena in the nerve membrane", Bull. Math. Biophys., Vol. 17, pp. 257-278, 1955. 49. R. Fitzhugh, "Thresholds and plateaus in the Hodgkin-Huxley nerve equa-tions", J. Gen. Physiol., Vol. 43, pp. 867-696, 1960. 50. R. FitzHugh, "Anodal excitation in the Hodgkin-Huxley nerve model," Biophys J., vol. 16, pp. 209-226, 1976. 51. K.S. Cole, R. Guttman, F. Bexanilla, "Nerve membrane excitation without threshold", Proc. Nat. Acad. Sci., Vol. 65, No. 4, pp. 884-981, 1970. 52. W.E. Boyce and R.C. DePrima, Elementary Differential Equations &  Boundary Value Problems, 3rd Ed., John Wiley & Sons, New York, 1977. 53. B. Frankenhaeuser, "Computed action potential in nerve from Xenopus Laevis," J. Physiol., vol. 180, pp. 780-787, 1965. 132 54. R. Fitzhugh, H.A. Antosiewicz, "Automatic computation of nerve excita-tion - detailed corrections and additions", J. Soc. Indust. Appl. Math., Vol. 7, No. 4, pp. 447-458, 1959. 55. V. Rowland, W.J. Maclntyre, T.G. Bidder, "The production of brain lesions with ele c t r i c a l current II - Bidirectional Currents", J. Neurosurg., Vol. 17, pp. 55-69, 1960. 56. H. Schimmel, H.G. Vaughan, "Design of a prosthetic system", in "Visual Prosthesis - The Interdisciplinary Dialogue", Proc of the 2nd Conf. on Visual Prostheses; Sterling, Bering, Pollack, Baughan, Eds., Academic Press (1971), pp. 236-239. 57. L. Lapicque, "Recherces quantitative sur 1*excitation electriques des nerfs traitee comme une polarization", J . Physiol., Paris, v o l . 187, pp. 622-635, 1907. 58. D. Noble, R.B. Stein, "The threshold conditions for i n i t i a t i o n of action potentials by excitable c e l l s " , J. Physiol., Vol. 187, pp. 129-182, 1966. 59. A.L. Huxley, W.A.H. Rushton, "The electrical constants of crustacean nerve fibre", Proc. R. Soc. B., Vol. 133, pp. 444-479, 1946. 60. R. Plonsey, Bioelectric Phenomena, McGraw-Hill (1969), pp. 19-20. 61. F. Offner, "Stimulation with minimum power", J. Neurophysiol., v o l . 9, pp. 387-390, 1946. 62. R. Plonsey, Bioelectric Phenomena, McGraw-Hill (1967), pp 80, Table 3.3. 63. Z.L. Kovacs, T.L. Johnson and D.S. Sax, "Estimation of the distribution of conduction velocities i n peripheral nerves", Comput. B i o l . Med., Vol. 9, pp. 281-293, 1979. 64. F. Buchthal, A. Rosenfalck and F. Behse, "Sensory potentials of normal and diseased nerves", in Peripheral Neuropathy, Vol. 1, P.J. Dyck, P.K. Thomas and E.A. Lambert, Eds., W.B. Saunders Co., pp. 442-464, 1975. 65. W.F. Collins, F.E. Nulsen and C.T. Randt, "Relation of peripheral nerve fibre size and sensation in man", Arch. Neurol., Vol. 3, p. 381, 1960. 66. B. Frankenhaeuser and A.L. Hodgkin, "The after-effects of impulses in the giant nerve fibres of l o l i g o " , J. Physiol., Vol. 131, pp. 341-376, 1956. 67. W.J. Adelman and Y. Palta, "The influence of external potassium on the inactivation of sodium currents i n the giant axon of the squid, loligo pealeli", J. Gen. Physiol., Vol. 53, pp. 685-703, 1969. 133 68. N. Moran, Y. P a l t l , E. Levitan and R. Stampfli, "Potassium ion accumula-tion at the external surface of the nodal membrane in frog myelinated fibres", Biophys. J., Vol. 32, pp. 939-954, 1980. 69. Y. P a l t i , N. Moran and R. Stampfli, "Potassium currents and conductance - comparison between motor and sensory myelinated fibres", Biophys. J., Vol. 32, pp. 955-966, 1980. 70. Y. P a l t i , R. Stampfli, A. Bretag and W. Nonner, "Potassium ion accumulation at the external surface of myelinated nerve fibres", Isr. J. Med. Sci., Vol. 9, p. 680, 1973. 71. E. Neher and H.D. Lux, "Rapid changes of potassium concentration at the outer surface of exposed single neurons during membrane current flow", J. Gen. Physiol, Vol. 61, pp. 385-399, 1973. 72. W.J. Adelman and Y. P a l t i , "The role of periaxonal and perineural spaces in modifying ionic flow across neural membranes", Curr. Top. Membranes Transp., Vol. 3, pp. 199-235, 1972. 73. W.H. Adelman and R. Fitzhugh, "Solutions of the Hodgkin-Huxley equations modified for potassium accumulation in a periaxonal space", Fed. P r o c , Vol. 34, pp. 1322-1329, 1975. 134 APPENDIX I THE FRANKENHAEUSER-HUXLEY EQUATIONS The model equations as reported by Frankenhaeuser and Huxley [1], and a l l model constants are repeated for convenience in this Appendix. The equation for membrane voltage: 51. - A fi - i - i - i - i 1 dt C L stim Na K P %) Ionic currents: (mA/cm2) . • ^ F 2 f [Na+]Q - [Na*], . ^ = PNa m h ^V + Vr ) Rf < (V-V JF/RT 1 - e r - 2 r , F2 , - V\ )F/RT IR = PK n2 (V+Vr) ^  ( (V-V JF/RT 1 - e r _ ^ p 2 r [Na+]o - e ( ™ r ^T Ip = Pp p2 (V+Vr) ^  ( (v-f-V JF/RT 1 - e r 135 The equations for m, n, p and h d T - .V1-*) - h1 The rate constants (ms-*) «h - -0.1 CV+10) / (i-e(v + 1 0>/ 6) Ph - 4.5 / (l+e<4 5-V>/ 1 0) affl - 0.36 (V-22) / ( l - e( 2 2~V ) / 3) P m - 0.4 (13-V) / (l-e<V-1 3>/ 2°) an = 0.02 (V-35) / ( l - e( 3 5"V ) / 1 0) 8n = 0.05 (10-V) / (1_e(V-10)/10) ap - 0.06 (V-40) / (i-e<4°-V ) / l°) Bp - -0.09 (V+25) / ( l - e( V + 2 5 ) / 2°) Constants: 136 P„ = 8*10~3 cm/s Na Fv = 1.2*10-3 cm/s K P = 0.54*10"3 cm/s P sodium permeability potassium permeability unspecific permeability g^ = 0.0303 mho/cm2 leakage conductance = 0.026 mV leakage current equilibrium potential [Na+]Q = 114.5 mM [ N a ^ = 13.7 mM [K+] = 2.5 mM = 120.0 mM external sodium concentration internal sodium concentration external potassium concentration external potassium concentration F = 96514.0 C/g/M Faraday's constant R = 8.3144 J/K°/M gas constant T = 295.18 Ke absolute temperature I n i t i a l conditions: m(0) = 0.0005 h(0) = 0.8249 n(0) - 0.0268 p(0) = 0.0049 V(0) = 0 dV dt t=0, I «, =0 ' stim 

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0096542/manifest

Comment

Related Items