Critical behavior in the Artificial Axon

The Artificial Axon is a unique synthetic system, based on biomolecular components, which supports action potentials. Here we examine, experimentally and theoretically, the properties of the threshold for firing in this system. As in real neurons, this threshold corresponds to the critical point of a saddle-node bifurcation. We measure the delay time for firing as a function of the distance to threshold, recovering the expected scaling exponent of $- 1/2$. We introduce a minimal model of the Morris-Lecar type, validate it on the experiments, and use it to extend analytical results obtained in the limit of"fast"ion channel dynamics. In particular, we discuss the dependence of the firing threshold on the number of channels. The Artificial Axon is a simplified system, an Ur-neuron, relying on only one ion channel species for functioning. Nonetheless, universal properties such as the action potential behavior near threshold are the same as in real neurons. Thus we may think of the Artificial Axon as a cell-free breadboard for electrophysiology research.


I. INTRODUCTION
Action potentials are an interesting dynamical system, but not easy to study due to the complexity of the neuron. We recently introduced the idea of producing action potentials in vitro; our cell free system is based on the reconstituted biological components (phospholipid membrane, voltage gated ion channels) and thus on the same microscopic mechanism for generating voltage spikes as real neurons; we call it the Artificial Axon (AA) [1,2]. While our ultimate goal is to study the behavior of simple networks of Artificial Axons (an initial demonstration is given in [3]), it is also interesting to develop the AA as a cell free platform for electrophysiology studies. Specifically, our hope is that AA platforms may one day replace, in the laboratory, cultured neurons, which typically come from rats. In addition, this experimental system provides motivation and inspiration for developing models, numerical experiments, and eventually theory in the general field of complex systems. In this context, we concentrate here on one feature of an individual AA, which has the allure of universality (and thus is relevant for real neurons too), but which is difficult to measure in real neurons, namely the behavior near threshold for firing an action potential. In the language of dynamical systems, the threshold for firing is related to the critical point of a saddle node bifurcation. This correspondence has of course been explored before in neurons [4][5][6][7][8][9]. Notwithstanding, the AA is a simpler system compared to a real neuron, both in practice and conceptually. It is a dressed down experimental system where we try to distill the minimal components which can still generate an action potential. Notably, it has only one ion channel species and relies on one ionic gradient across the membrane which is the active element of the system. The role of a second ionic species and channel type is played instead by a simple electronic device which we call a current limited voltage clamp (CLVC) [1]. The command voltage to this clamp * zocchi@physics.ucla.edu is a control parameter for the dynamical system slightly different from the control parameters in traditional electrophysiology. As a function of this control parameter the system exhibits critical slowing down and scaling of the time to fire. Our purpose with this paper is to use the AA as a conceptual tool -a model -to generate the simplest interesting dynamical systems related to action potentials. We start from some experimental measurements of critical behavior in the physical AA and then move on to discuss this critical behavior for some corresponding dynamical systems. But first, let us introduce the AA, which we have described before [1,2]. It is a black lipid bilayer about 200 µm in size, separating two compartments ("in" and "out") with ionic solutions at different salt concentrations (and matched osmotic pressure), see Fig. 1.
Typically we keep the inside chamber at a salt concentration [KCL] = 30 mM and the outside at [KCL] = 150 mM . Given a small leak conductance for K + ions across the membrane (which is impermeable to other ionic species), equilibrium with respect to ionic currents is reached for a positive Nernst potential inside, with respect to the grounded outside: where |e| is the charge of the K + ion, T the temperature, k the Boltzmann constant, and square brackets denote concentration. Embedded in the phospholipid membrane are order of 100 voltage gated potassium ion channels: these are the active elements of the system. Voltage gated ion channels are the ionics equivalent of transistors in electronics; they are membrane proteins which act as pores for one specific ion species (potassium in our case), are impermeable to other species, and switch from closed (not conducting) to open (conducting) depending on the voltage across the membrane. We use an Archaean channel (KvAP) which can be expressed in E. coli and reconstituted in vesicles [10], which are then fused to the AA membrane (see Mat. & Met. voltages below about −50 mV and opens at positive voltages; its main characteristics are that it is slow (opening rates at positive voltage are of order ∼ 25 ms [1]) and that at positive voltage it eventually enters an "inactive" state which is functionally closed, but from which the open state is not accessible. Recovery from this inactive state necessitates negative voltages (< −100 mV ) and is also slow (measured in the hundreds of ms). The detailed channel dynamics are complex [10], but we will not need to enter into all the details for this study. Fig. 2 shows the measured open probability curve p o (V ) for the channels in our system, reproduced from [11]. The meaning of this plot is that it represents the equilibrium probability that the channel is open if the voltage was held steady at V , and if the inactive state was not there. The plot is obtained by quickly stepping the voltage to V (using a voltage clamp) from a large negative value where all the channels are in the closed state, and measuring the peak current (the current through the channels before they inactivate). The last component of the AA is the current limited voltage clamp (CLVC). This simple device consists of a traditional voltage clamp with a large resistance in series ( Fig. 1). Given a voltage V (t) in the AA (in the inside chamber, relative to the grounded outside chamber), a command voltage V clamp to the CLVC causes a current role of the CLVC is to hold the system off equilibrium at a negative "resting potential" V r , while still allowing the voltage in the AA to fluctuate in response to perturbations. We can now write down the equation governing the voltage in the AA (basically the same equation is the starting point for all models of the neuron [9]): (2) C is the membrane capacitance, N 0 the number of ion channels, χ the single channel conductance (with the channel open), χ a small leak conductance which is present even if the channel is closed (χ << χ), p o (t) the probability that the channel is open, V N the Nernst potential (1), R the series resistance of the CLVC, and V clamp the command voltage to the CLVC. If we multiply both sides by C, the equation describes the charging of the membrane capacitance by the currents specified on the RHS: the current through the channels, with driving force proportional to [V N − V (t)], and the current injected by the CLVC, with driving force proportional to V clamp −V (t). In general the probability p o (t) is specified by a set of rate equations which reflect the dynamics of the channels, a description introduced by Hodgkin and Huxley in their famous 1952 paper [12]. We come back to this point in Section III. For now we note that if we hold V clamp at a large negative value, a steady state solution of (2) exists such that V (t) = V r is also large and negative, and therefore p o ≈ 0 (see Fig. 2), namely: From this expression we see that the CLVC fills the role of a second ionic gradient in determining the resting potential, with V clamp corresponding to the Nernst potential for this second species (e.g. Na + ) and 1/R corresponding to the leak conductance times the number of channels. The CLVC resistance R is chosen such that the clamp current is sufficient to pull the resting potential to negative values V r ∼ −100 mV with channels closed, and on the other hand, such that the channel current with channels open can overwhelm the clamp current so that the AA can fire. From (2) and (3) one finds that we require N 0 χ < 1/R << N 0 χ. Typical values of the parameters for the AA with the KvAP channels are: C ≈ 300 pF , N 0 ∼ 100, χ ≈ 10 pA 50 mV = 2 × 10 −10 Ω −1 , As a result, V r ≈ −90 mV , currents are in the 100 pA − 1 nA range, and the characteristic time scale C/(N 0 χ) ∼ 10 ms while C/(N 0 χ ) ∼ 1 s. Under these conditions, starting from the resting potential V (t) = V r , a perturbation such as a positive input current (which, in the experiment, can be delivered through a separate current clamp) can cause the system to fire an action potential [2]. Referring to Fig. 3, the rising edge is caused by opening of the potassium channels, and has a universal shape; the falling tail is caused by the channels going into the inactive state, which allows the CLVC to pull the voltage down. Elsewhere we will discuss the requirements on the channel dynamics in order for the system to exhibit self-sustained oscillations (that is, a firing rate) in the presence of a steady input current; here we focus on the threshold for firing, which is a critical point of the corresponding dynamical system.

II. EXPERIMENT
We have seen in Section I that the AA is a much simplified version of a spiking neuron. Nonetheless, there are aspects of the dynamics of the neuron which possess a degree of universality, so we expect them to remain the same independent of the details of the system. In particular, the behavior near a critical point, i.e. the threshold for firing. To explore the neighborhood of this critical point, we use as control parameter the command voltage to the CLVC (V clamp ). Physically this differs slightly from traditional electrophysiology, as the CLVC is neither a voltage clamp nor a current clamp, but rather something in between. In a neuron where the action potential is shaped by sodium and potassium currents, the analogous control parameter would be the Nernst potential of the potassium ions. The experimental protocol is as follows. We start from a steady state with V clamp held at −200 mV , which typically corresponds to V r ≈ −180 mV . At these voltages the channels are fully closed, but not inactivated. The clamp command voltage V clamp is then stepped up to various values in the range between 0 to −100 mV , and held steady. The consequent increase in the AA voltage V (t) is analogous to charging an RC circuit, with the membrane acting as the capacitor and the leak current and clamp resistor acting as the resistances.
If functional ion channels are in the membrane, and V clamp is above a threshold V crit , the AA will fire, i.e. a voltage spike will occur. The current through the membrane due to the channels opening is much larger than the CLVC can provide, which causes the voltage to spike towards the Nernst potential. If V clamp is below V crit , the AA does not fire and voltage will stabilize at V (t) ≈ V clamp . After the AA has fired, V clamp is stepped back to −200 mV and held there for ∼ 20s to allow the channels to recover from inactivation; then the system is probed again through another cycle. One challenge in these experiments is the stability in terms of maintaining the same threshold over a set of measurements meant to explore the vicinity of the critical point. Thresholds may change over time due to changes in the number of channels in the lipid bilayer, and changes in the capacitance (due to expansion or shrinking of the lipid bilayer).

A. Firing Time Delay
When V clamp is stepped above threshold, the KvAP channels will open to initiate a spike and then inactivate, preventing further spikes until they are allowed to recover at resting potentials < −120 mV . The spikes that occur above threshold are qualitatively similar regardless of V clamp , a key feature that defines action potentials. Instead, the value of V clamp determines the amount of time that elapses before the onset of the action potential. Namely, we find that as the membrane voltage is stepped closer to the voltage threshold V crit from above, there is a significant delay in the time it takes to depolarize the inner chamber. The measure of delay τ we adopt here is the time interval between the step in V clamp and the peak of the action potential. This differs slightly from the time of spike onset, with the main benefit being the ease of locating the peak versus the threshold. In any case, both methods produce equivalent results in what follows. We have also explored alternative definitions, such as measuring the first derivative peak, with similar results for the scaling of τ .

B. Determining the Threshold
We expect the delay time τ to scale with the distance to the critical point (V clamp − V crit ), as explained in Section III. To experimentally confirm this relation, we must first locate V crit with precision. Unfortunately, there is no uniform method to locate this threshold and it is generally accepted to be an empirical quantity [13], though some models do attempt to calculate it numerically [14]. We take a slightly different approach given the fact that V crit is the critical point of our system. To find V crit , we assume the existence of this critical point and proceed as follows. We take a series of action potentials obtained from the same AA prep and measure the delay time τ . The delay is subsequently plotted against (V clamp −V crit ) on a log-log plot, with an initial estimate for V crit being the empirically observed value during experiment. We then adjust V crit to make the points fall on a straight line. This is a usual procedure in the study of critical phenomena, and we find that it gives a more precise -and different -determination of the critical point compared to other methods used in electrophysiology, which we summarize later. In detail, we produce, for the same data set, log-log plots of τ vs (V clamp − V crit ), for different values of V crit . We fit each plot to a straight line, subtract the value of the fit from the experimental data, and fit a quadratic form to these residuals. We then determine the value of V crit which corresponds to the coefficient of the quadratic term vanishing (Fig. 4).
If we want to measure critical exponents, a precise determination of the threshold is essential, where by threshold we mean the location of the critical point in parameter space. The methods that have been employed in electrophysiology for estimating action potential thresholds were also considered. The well known methods all relate to observing changes in the derivative of the voltage [15]. From our perspective, these methods have a tendency to overestimate the voltage threshold, due to the fact that channels have a slight delay in opening and by the time there is a significant peak in the derivative the voltage has already risen beyond the threshold. Indeed, the location of the threshold shown in Fig. 3 is not where one would expect from a qualitative examination of the time traces. It does however correspond to the location of the critical point in parameter space, as we show in Section III. For similar recordings and threshold determination in real neurons, see for instance [16,17]. The result of the scaling assumption for one data set is shown in Fig. 5. The slope of the fitted straight line is −0.53, close to −1/2. Unfortunately the range of the measurements is quite limited, covering only a decade in the control parameter (V clamp − V crit ). As we will see in Section III, the relatively small number of channels in the AA limits how close one is able to get to the critical point, and beyond that, the stability and noise in the system. Fig. 6 shows plots for two more data sets (independent preps of the AA). These data sets show differing values of the scaling exponent, we suspect due to the drifting of the threshold during the experiment as mentioned earlier. In particular, excess channels that do not initially insert into the membrane may slowly insert over time. This causes the number of channels in the membrane to increase over time which alters the threshold (see section III and Mat. & Met.).

C. Materials and Methods
The artificial axon (AA) system consists of voltage gated ion channels from the Archaea Aeropyrum Pernix (KvAP) inserted into a phospholipid bilayer (black lipid 3-phosphocholine dissolved in chloroform at 25 mg/mL (Avanti Polar Lipids) is Nitrogen dried and redissolved in 12.5 µL decane for a concentration of 20 mg/mL. A small droplet of the decane-lipid mixture is dripped via syringe onto the aperture of the aforementioned plastic cup, which serves as the intracellular chamber (Fig. 1). The cup is placed under vacuum for 30 minutes to ensure that the chloroform is thoroughly removed before the cup is inserted into the Teflon chamber. The outer chamber is filled with 150 mM KCl and 10 mM Hepes pH 7.0 while the plastic cup is filled with 30 mM KCl and 10 mM HEPES plus 120 mM sucrose pH 7.0 to equalize osmotic pressure. Additional DPhPC is pipetted in at the bottom of the cup to seal the opening and prevent the diffusion of ions. At this point the electronics are turned on and the CLVC injects current into the inner cup which clamps the voltage at the desired value. All experiments start with the voltage clamped at -200 mV after the system is prepared. A planar lipid bilayer is then painted into the opening and allowed to stabilize. Expressed KvAP solubilized in DPhPC vesicles are then removed from an −80 • C freezer and thawed quickly at room temperature. This solution is mixed with a pipette prior to a small amount (.3-1.5 µL) being injected close to the lipid bilayer. The channels are given a few minutes to insert, then the CLVC is stepped up from its initial value of -200 mV to confirm the insertion of channels. The occurrence of an action potential signifies that a sufficient number of channels have inserted into the membrane. A typical voltage trace from which the measurements are obtained is displayed in Fig. 7. The ion channels themselves were expressed in E coli. following published procedures [1], which were adapted from [18]. The procedure for reconstitution of the channels into lipid vesicles also follow that of [1], with the adjustment of the lipid to protein ratio to 3 instead of 1 (wt/wt).

III. THEORY
Models of the neuron as a dynamical system [9] hail back to the ground breaking studies of Hodgkin and Huxley [12] and form some of the most quantitative descriptions in neurobiology. The AA is a much simplified version of a neuron, a minimal system on which some of the universal characteristics of the more complex living matter can be illustrated and studied. Such is the threshold for firing an action potential, which is one of the most important and defining characteristics of neurons. It is not difficult to see that, in the AA, the threshold for firing corresponds to the critical point of a saddle node bifurcation. We show this in subsection B, after briefly recalling the general theory [19] in subsection A. Then we discuss critical behavior in the simplest model with channel dynamics similar to the AA (subsection C).

A. Normal form of the saddle-node bifurcation
Consider the 1D dynamical system: where r is a control parameter. From the graph ofẋ vs x (a parabola) we see that for r < 0 the system has two fixed points (ẋ = 0), one stable and one unstable. They merge at the critical point r = 0; for r > 0 the velocity is always positive so the system escapes to infinity. However, near the critical point it exhibits critical slowing down: for r > 0 (4) integrates to: We see from (5) that the time to escape to infinity is finite; starting from a large negative value x(0), so that b ≈ −π/2, that time is: which can also be estimated by: The delay time diverges as one approaches the critical point from above, and scales with the distance to the critical point with characteristic exponent −1/2.

B. Saddle-node bifurcation in the Artificial Axon
The arguments in this section closely follow those developed in [4]. Consider the "membrane equation" (2) for the AA. To show the existence of the bifurcation, we make the approximation of considering the channel dynamics fast enough that the probability that channels are open p o (t) in (2) always is the equilibrium value p o (V ) shown in Fig. 2, and we ignore the inactive state. The exact meaning of this approximation is seen more clearly in the next section. In this case, (2) becomes the 1D dynamical system: with given values of the parameters q and V 0 (kT is the thermal energy), as shown in Fig. 2. The fixed points (dV /dt = 0) of the dynamical system (8), (9) are given by: Fig. 8 shows plots of the LHS and the RHS of (10), the relevant range of voltages being V clamp < V (t) < V N . We see that for V clamp sufficiently negative (Fig. 8a) there are 3 fixed points. As V clamp is increased, the stable and unstable fixed points on the left merge at a critical point (Fig. 8b) where the two curves in the figure are tangent to each other. For larger values of V clamp these two fixed points disappear and only the stable fixed point close to the Nernst potential remains. Evidently this is the same bifurcation as exhibited by the dynamical system (4), so we expect the same behavior near the critical point. The disappearance of the lower stable fixed point for V clamp > V crit corresponds to firing of the action potential. If V clamp is stepped close to, but above, V crit , there will be a delay time for firing, scaling as: Explicitly, the critical point Set V clamp = V crit + and write (8) in the form: Expanding around the critical point: Using (10) one recognizes that the first two terms on the RHS of (15) vanish. Also, the coefficient of the quadratic term is positive. Finally, Thus, close to the critical point the dynamical system (8) reduces to (4). Using the estimate (7) for the delay time: and, noting that b ∝ N 0 χ/C: The delay time scales with the distance to the critical point with the exponent −1/2, with a prefactor proportional to the characteristic RC time scale, modulated by a factor which depends on the number of channels and the ratio between the open channel conductance and the CLVC conductance. The rate at which the voltage proceeds to the fixed point close to the Nernst potential can be visualized in a phase space plot. Fig. 9 shows a plot of dV /dt vs V for several trajectories (time traces) obtained from integrating the AA model described in the next section. When V clamp is stepped up, dV /dt jumps to a positive value and the voltage immediately starts to rise. As V approaches the critical point V 1 , dV /dt will reach its minimum value, which can be arbitrarily small depending on how close V clamp is to V crit . The scaling of the delay time τ can be attributed to this bottleneck region. Eventually, dV /dt will start to rise again and spiking will occur, as shown in the right half of Fig. 9. The subsequent drop and sign change in velocity is due to the inactivation of the channels which moves the fixed point back near V crit . The theory in this section is not new; specifically, this dynamical systems approach appears in the excellent article by Morris and Lecar [4], where they describe their electrophysiology measurements on the barnacle muscle fiber. We have presented the calculations in detail to show the generality as it pertains to the AA. In summary: for the delay time τ measured in the experiments of Section II, we expect a scaling exponent of −1/2. However, (18) also shows that in order to see that scaling, the experiment must be stable with respect to the parameters appearing in the prefactor to the power law. If, during the course of one run, the capacitance C changes (because the membrane shrinks), or the number of channels (and thus the parameter b) changes (for example, some channels may diffuse to the rim of the membrane and effectively disappear; new channels may insert), the exponent measured from the log-log plot of the data may be affected. We believe such a drift affected the data sets in Fig. 6, which show an exponent significantly different from −1/2 (discussed in section IV).

C. Model with channel dynamics
In order to represent the complete action potential seen in Fig. 3 -not just its rising front -we need to introduce channel dynamics in the AA model, specifically the inactive state. The minimal model for the KvAP channel has 3 states: closed, open, and inactive (which is functionally closed), and voltage dependent rates to move between these states. The actual dynamics of the channels -if one looks in detail -is significantly more complex [10], with more states and rates. For instance, the KvAP is a tetramer, and the conformational state of each subunit must in principle be specified for a more complete description of the state of the channel. However, these microscopic details are often not important for the description of the kind of collective phenomena -such as action potentials -that we want to explore with the AA. In addition, a more complete description leads to the introduction of several more rates, i.e. an uncomfortable proliferation of the number of parameters in the model. We therefore start from a minimal model of the chan-nels which still contains the important dynamics: it is depicted schematically in Fig. 10(a) Fig. 10b. It is not immediately clear which of these two models is a better representation for the KvAP; for our discussion of scaling near the critical point they lead to the same result for universal quantities such as the scaling exponent, but other important behavior may differ. We also note that if the time scales 1/k o , 1/k c are much faster than any other time scale in the system (this is not the case for the present AA), then the dynamics of the two models are equivalent. In the following, we will stay with Model A. Since the channels are voltage gated, the rates k o , k c must be made voltage dependent. In the usual way, we assume a 1D barrier crossing process and write: This "symmetric" form contains the extra assumption that the transition state lies exactly "in the middle" between the open and closed states along the 1D conformational landscape; it can be relaxed at the cost of introducing one more parameter. The rates (19)  We use values for the parameters which are based off of the experimental data, by fitting the model below onto the voltage traces in the experiment; they are: The AA model which includes channel dynamics (according to Model A) consists of equation (2), which we re-write here: (21) and the following rate equations for the probability that the channels are in the open state (p o (t)), closed state (p c (t)), or inactive state (p i (t); see Fig. 10): which enforce p o + p c + p i = 1; the two equations (22) are coupled to (21) through the voltage dependence of the rates (19), (20). The 3D dynamical system (21), (22) reproduces the phenomenology of the physical AA. Fig. 11, which is analogous to Fig. 3, shows two action potentials evoked in the model by stepping V clamp to different values. We used this model to simulate the experiments of Section II, i.e. we produced a series of action potentials for different values of V clamp , measured the delay time τ , and determined the threshold V crit using the method described in Section II. Fig. 12 shows the resulting log-log plot of τ vs (V clamp − V crit ). We used the following parameter values in the model: C = 329.7 pF, χ = 1.67 × 10 −10 Ω −1 , χ /χ = 8.8 × 10 −4 , R = 2 × 10 9 Ω, N 0 = 110 , and the initial condition V (t = 0) = −200 mV and p o (t = 0) = p i (t = 0) = 0 (i.e. channels closed initially). The advantage of the simulation is of course that we obtain a larger range of the power law behavior and thus a good determination of the scaling exponent; the linear fit in Fig. 12 has a slope of −0.51. Here there is no question of the stability of the "experiment", since parameters such as C and N 0 are fixed, but what limits how close one can get to the critical point (and thus how large a delay time τ one can observe) is the finite number of channels.

IV. DISCUSSION
Simulation of the model in section III allows us to explore regions close to the critical point that are difficult to access in the experiment. Fig. 12 shows that the simulation agrees well with experiment in the range where they overlap. However, one limit still in effect is the finite size of the system, namely the number of channels. Given that the number of channels in the simulation must match the experiment from which the rates were taken, 13. Numerically calculated relation between the number of channels in the system and the threshold, using steady state opening probabilities and no leak. The blue fit is logarithmic with coefficient −kT /q (see eq. (9) and [13]). the scaling exponent will only be correct up to a certain resolution, and we cannot achieve arbitrarily long delay times. The simulations also demonstrate the 'stability' of the scaling exponent, in that significantly varying the rates and conditions do not alter the scaling behavior correspondingly. This is of course to be expected since the scaling phenomenon arises from the saddle-node nature of the critical point. Additional simulations of model A were performed to determine how the number of channels and the leak conductance affect the threshold. Solutions to eq. (10) with steady state opening were also calculated numerically to obtain a relation between N 0 and V crit for varying values of χ . The results show that there exists a range of channels number N 0 such that excitability is possible, above and below which there is no bifurcation and thus no possibility of action potentials occurring. This fact can also be seen directly from the graphical representation of equation (10) in Fig. 8 : increasing N 0 makes the slope of the RHS of (10) shallower, and thus V crit more negative, but since χ /χ > 0 there is a maximum N 0 beyond which the only fixed point is the one close to V N (i.e. a steady state with channels open). Conversely, decreasing N 0 makes the slope of the RHS of (10) steeper, and thus V crit increases; the critical point V 1 (eq. (12)) moves towards the inflection point of the open probability curve p o (V ), and we see that if N 0 is too small, the unstable fixed point disappears, and so there is no firing. In the case that there is no leak (χ = 0), we found a logarithmic relationship between the number of channels and the threshold. This is in agreement with previous numerical studies [13] calculating a relationship between the conductivity and the threshold voltage, since the total conductance is directly proportional to the number of channels. For larger values of the leak, we find that V crit is significantly more negative and the relationship is no longer logarithmic. When the leak is taken to be the same value as the simulation, the numerical solutions are in good agreement with the simulation of model A, indicating that the threshold does not depend strongly on the opening model. These results are displayed in Figs. 13 and 14. It is interesting to note that if the data from Fig. 6 are taken and plotted with the simulated threshold rather than the fitted threshold, a scaling exponent much closer to the theoretical one is obtained. This further supports our reasoning for the deviation in experimental behavior (i.e. the drift in voltage threshold as the experiment progressed). We expect that if a drift had occurred, it would have decreased the voltage threshold over time, as channels floating above the lipid membrane would have inserted over time, increasing the number of channels in the system. As such, a simulation using fits on the first few data points would be well suited to provide an estimate for how the system would behave in the absence of this drift. Our simulation rates are fitted onto the spikes closest to the critical point, which were the data points first recorded in the experiment. In conclusion, we report measurements of critical behavior in the AA, and corresponding simulations of a minimal model which describes the system well. Critical properties are universal and are thus the same for the AA, the model, and real neurons. Our quest is to construct, at each scale, the minimal system which retains the essential properties of the more complex biological systems of today. In this (primordial) world, the AA is the minimal neuron, networks of artificial axons may be developed into minimal nervous systems. More prosaically, the present work is an example of how the AA can be developed into a cell-free breadboard for electrophysiology research.