The potential importance of saturating and supersaturating contrast response functions in visual cortex

Abstract

Most cortical visual neurons do not respond linearly with contrast. Generally, they show saturated responses to stimuli of high contrast, a feature often characterized by a divisive normalization function. This nonlinearity is generally thought to be useful in focusing the dynamic response range of the neuron on a particular region of contrast space, optimizing contrast gain. Some neurons not only saturate but also supersaturate; at high contrast, the response of the neuron decreases rather than plateaus. Under the contrast gain control theory, these cells would seem to reflect a nonoptimal normalization pool that provides excessive inhibition to the neurons. Since very few data on supersaturation are available, this article examines the frequency with which such neurons occur in macaque visual cortex by considering an extension of the Naka–Rushton equation with the capacity to represent nonmonotonic functions. The prevalence of gain-control theories for saturation has occluded an additional computational function for saturation, namely, in detecting the conjunction of certain features. A saturating nonlinearity is a critical part of the selective detection of compound stimuli over their components. In this role, the existence of saturating contrast response functions might be considered necessary rather than simply optimal.

Introduction

Most cortical visual neurons in cat area 17 and macaque areas V1 and V2 (Albrecht & Hamilton, 1982; Maffei & Fiorentini, 1973)do not respond linearly with respect to the contrast of the stimulus on their receptive fields. A few cells show some degree of expansive nonlinearity at low contrast, which might serve to sharpen tuning (Heeger, 1992). Most show a compressive nonlinearity at high contrast such that their response saturates or plateaus. 

Although saturation might represent some basic limitation on the neuron's biophysical capabilities to respond more strongly, this seems unlikely for several reasons. Whereas some cells might plateau at a response rate of 30 impulses per second (ips), other neurons are fully capable of responding well in excess of 200 ips. In fact, even a single neuron will saturate at different response rates, dependent on other characteristics of the stimulus. For example, a cell that saturates at some given contrast level for optimal stimulus orientation will, for non-optimally oriented stimuli, typically saturate at the same contrast but reach a lower peak response rate (Albrecht & Hamilton, 1982). Therefore, contrast response saturation seems to result from a process that actively controls the response rate of a cell. This might be achieved by a set of inhibitory connections from a pool of neurons, the normalization pool, that increasingly inhibit the neuron as contrast rises (Carandini, Heeger, & Movshon, 1997; Heeger, 1992). The normalization process results in saturating contrast response functions. This may effectively control the gain of the contrast response function, and may serve to reduce redundancy in the natural image (Schwartz & Simoncelli, 2001). 

There is evidence that optimizing the dynamic range of the neuron is at least one of the functions of saturation. Neurons alter their contrast response functions according to recent levels of contrast (Ohzawa, Sclar, & Freeman, 1985), enabling them to center the rising slope of the contrast response function on the appropriate range. Also, the contrast response functions across the population of fly visual neurons seem to be optimized for the range of contrasts found in natural scenes (Laughlin, 1981), as are the functions in cat LGN and the magnocellular neurons of macaque LGN (Tadmor & Tolhurst, 2000). In macaque V1, however, this is not the case (Clatworthy, Chirimuuta, Lauritzen, & Tolhurst, 2003). Here, neurons respond optimally to a range of contrasts higher than generally observed in the natural environment and show highly variable degrees of saturation. 

An additional problem for the hypothesis that visual cortical neurons saturate to optimize their dynamic range is that some V1 cells “supersaturate”. That is, their responses actually decrease at very high stimulus contrasts. This behavior has been noted by numerous authors (Albrecht & Hamilton, 1982; Bonds, 1991; Ledgeway, Zhan, Johnson, Song, & Baker, 2005; Li & Creutzfeldt, 1984; Maffei & Fiorentini, 1973; Mancilla, Fowler, & Ulinski, 1998; Tyler & Apkarian, 1985) and has been simulated in a circuit model of V1 neurons (Somers et al., 1998). However, it has not been properly quantified, and its computational implications are rarely considered (but see Ledgeway et al., 2005). One reason that supersaturation has received relatively little attention is our tendency to characterize response functions with the equation used by Naka and Rushton (1966), which is incapable of representing nonmonotonic data. 

The aims of this article are twofold: to characterize the contrast response curves of neurons in V1 and V2 of the macaque using an adaptation of Naka and Rushton's equation, such that supersaturation can be appropriately represented, and to reconsider the role that the nonlinear responses of these neurons (both saturating and supersaturating) may have in visual processing. 

Experiment 1: Characterizing supersaturation

Methods

Collection of neuronal data

Recordings were taken from 541 neurons in the visual cortex of 31 anesthetized, paralyzed macaques as part of a larger series of experiments. Most of the data shown here represent a reanalysis of data that have been previously published (Solomon & Lennie, 2005; Solomon, Peirce, & Lennie, 2004; Webb, Dhruv, Solomon, Tailby, & Lennie, 2005) or are in preparation for publication. Details of the surgical preparation are described in greater detail elsewhere (Solomon & Lennie, 2005). All procedures conformed to the guidelines approved by the New York University Animal Welfare Committee. 

Anesthesia was induced with ketamine hydrochloride (Vetalar, 10 mg/kg ⁻¹) and maintained during surgery with thiopental sodium. The monkey was intubated, the head placed in a stereotaxic frame, and a craniotomy made over the occipital cortex, centered on or near the lunate sulcus. Postsurgical anesthesia was maintained by continuous infusion of sufentanil citrate (4–12 μg·kg ⁻¹·hr ⁻¹) in physiological solution (Normosol-R, Abbott Laboratories, Illinois, USA) with added dextrose (2.5%). Muscular paralysis was then induced and maintained by continuous infusion of vecuronium bromide (100 mg·kg ⁻¹·hr ⁻¹). The monkey was respirated artificially to keep end-tidal CO ₂ near 33 mm Hg. Electroencephalographic and electrocardiographic data were monitored continuously throughout surgery and the experiment. Rectal temperature was kept near 37°C using a heating blanket. Recordings were made using either tungsten-in-glass (Alan Ainsworth Electrodes, UK) or epoxy-coated tungsten electrodes (FHC Inc., Maine, USA) with an impedance of 1–5 MΩ. Spikes were isolated from single neurons using a combination of thresholding and template-matching techniques. 

Histological track reconstructions showed that 263 of the neurons were in V1 and 205 were in V2. For the remaining 73 cells, the area could not be precisely determined (e.g., because they fell too close to the V1/V2 border), and those cells are not considered further in this study. 

For each neuron, the optimal stimulation parameters were determined for circular patches of drifting sinusoidal grating varying in location, size, spatial frequency, and orientation. The contrast response function was then measured for each neuron using optimal stimulus parameters for the drifting gratings. This was generally conducted with an achromatic grating, unless the neuron responded poorly to achromatic stimuli, in which case a sinusoidal grating in the neuron's preferred color was used. The vast majority neurons were probed with a blank screen and six Michelson contrast levels ranging from 0.031 to 1.0 in 1-octave steps. For a few cells, additional contrasts were tested. Cells' responses were analyzed in terms of either the response rate in impulses per second ( f ₀) or the amplitude of modulation of that response at the frequency of the stimulus drift rate ( f ₁), whichever was larger. 

Modeling contrast response functions

Most groups (e.g., see Albrecht & Hamilton, 1982; Sclar, Maunsell, & Lennie, 1990) characterize contrast response functions using the hyperbolic ratio equation of Naka and Rushton (1966), which is closely related to the Michaelis–Menten equation for enzyme kinetics¹ (Michaelis & Menten, 1913): 

$R=R_{\max }\frac{c^n}{c_{50}^n+c^n}+B,$

(1)

where R is the output response of the neuron, B is the baseline response of the neuron, R_max represents the maximum response of the neuron, c₅₀ represents the contrast at which the response is halfway between baseline and maximum, nis simply referred to as the exponent, and c is the contrast presented to the neuron. 

To model the response of a neuron that is capable of supersaturating, we must add an additional parameter:  

$R=R_{\max }\frac{c^n}{c_{50}^{san}+c^{sn}}+B,$

(2)

where s is an additional parameter allowing the suppressive exponent to vary at a different rate to the excitatory exponent. The traditional Naka–Rushton equation is the specific form of this, with s = 1. The characteristics of the function are similar to Equation 1 in that the exponents control the shape of the function, R controls its amplitude, and c controls its contrast scaling. Unfortunately, with s ≠ 1, R _max and c ₅₀ values do not relate simply to the intuitively appealing maximal response and half-maximum contrast of the function, respectively. The R _max value remains the point at which the curves asymptote, although for nonmonotonic forms, this is obviously less than the maximal response. The c ₅₀ value for nonmonotonic functions represents the contrast at the curve's “shoulder” rather than at its half-maximal response. Examples can be seen in Figure 1. 

Figure 1

 

The form of the modified Naka–Rushton function as the suppressive exponent, s, progresses from 1 (equivalent to the original Naka–Rushton equation) to 1.7, with all other parameters fixed ( R _max = 100, c ₅₀ = 0.3, n = 2).

Data from the population of neurons were fit with both functions using Matlab's fmincon() function, aiming to minimize the chi-square fit term given by  

${\chi }^2=\sum _i\frac{(e_i-o_i{)}^2}{{\sigma }_i^2},$

(3)

where i is the index of this particular contrast level, e is the expected response at this contrast level given the current model parameters, o is the observed response, and σ ² is the trial-wise variance in responses at this contrast. Where there was zero variance in responses (arising from no response on any trial), the σ ² was set to 0.001 to prevent the error term from expanding to infinite. In fitting, the exponent terms of the models were bounded in the region 0.3–4.0. 

There is a problem in trying to test which of a pair of models with different numbers of free parameters fits a dataset “better”. One general solution is to measure how well the data are fit by the model and then impose some form of penalty on models for having more parameters. Aikake's Information Criterion (AIC; Motulsky & Christopoulos, 2004) is one commonly used method. The AIC value for a model fit to the data is calculated from 

$\mathrm{A}\mathrm{I}\mathrm{C}=n\log (\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{)}\mathrm{+}\mathrm{2}\mathrm{p}\mathrm{,}$

(4)

where n represents the number of data points being modeled, p represents the number of parameters in the model, and RMS is the root-mean-square error of the model fit. The optimal model under this method is taken to be the one with the lowest AIC value (nearest to negative infinity). A related method is the Bayesian Information Criterion, which is identical except that the penalty term is log(n)p instead of 2p. In this study, there were nearly always seven contrast levels, and because log(7) ≈ 2, the methods are equivalent. 

An alternative method, referred to as the normalized chi-square by Cavanaugh, Bair, and Movshon (2002), takes the χ² error term (Equation 3) and normalizes by the degrees of freedom for the model (Hoel, Sidney, & Stone, 1971): 

${\chi }_N^2=\frac{{\chi }^2}{df}.$

(5)

As with the AIC, the best fitting and most efficient model is taken to be that with the lowest χ_N². In practice, for the given data and models, this method typically inflicts a more stringent penalty on models with additional parameters. 

Quantifying supersaturation

Ledgeway et al. (2005) quantified the degree to which the response of neurons rises monotonically with respect to contrast using a monotonicity index (MI) defined as 

$\mathrm{M}\mathrm{I}=1-\frac{R_{\max }-R_{100}}{R_{\max }-R_0},$

(6)

where R_max is the maximum response of the neuron, R₁₀₀ is the response at maximal contrast, and R₀ is the baseline response. This index takes a value of 1 where the neuron's response rises monotonically and less than 1 for neurons demonstrating supersaturation. A neuron whose response at the maximal contrast falls fully back to baseline rate takes an MI of 0. In the current study, the MI was calculated for both the raw data and the model fit to the raw data. 

Results

Comparing the quality of fit for the models

Figure 2 shows a comparison of model fits for a series of sample neurons ranging in contrast response from very weak saturation to extreme supersaturation (the behavior shown in Figure 2f is unusual). For the neuron with weak saturation ( Figure 2a), the Naka–Rushton equation characterized the neuron's response function very well, such that the additional parameter of the modified version was unnecessary. This subjective decision is in agreement with the χ _N ² error term, which is lower for the Naka–Rushton equation (0.106) than for the modified version (0.117). The neurons shown in the remaining panels were all better fit by the modified Naka–Rushton model. In the supersaturating functions, it is clear that this is because of the inability of the Naka–Rushton function to provide a nonmonotonic response. However, the cell shown in Figure 2b shows a behavior that was also common, where the neuron did not show any sign of supersaturation but still showed a significantly improved fit with the modified model. This type of behavior can often only be shown for neurons with very consistent responses (note that the error bars in Figure 2b are almost obscured by the data points). The AIC error term, which is less stringent on models with additional parameters, was smaller for all of the six neurons shown in Figure 2. 

Figure 2

 

Comparing the fit quality of the standard Naka–Rushton function (broken line) with the modified form (solid line) for various cortical neuron samples with varying degrees of saturation.

For both models, the quality of the fits was extremely high across the population of neurons; the Naka–Rushton model explained, on average, 96.9% of the variance in responses and the modified model explained 98.7% of the variance. According to the AIC error term, the modified function performed better than the Naka–Rushton function in 175 of 264 V1 neurons (66.3%) and in 141 of 205 V2 neurons (68.8%). Even the normalized chi-square term suggests that it was an improvement for nearly half of the neurons tested (47.2% of V1 neurons and 45.9% of V2 neurons). 

Occurrence of supersaturation

Measuring its occurrence from the raw responses of the neurons, we find supersaturation in 45 of 263 V1 neurons (17.1%) and in 53 of 205 V2 neurons (25.9%). This probably represents an underestimate of the frequency with which it occurs due to the relatively sparse sampling at high contrasts. In particular, cells may well peak between 50% and 100% contrast, between which no stimuli were generally presented. Using the best model fit for the neuron, we find supersaturation in 65 neurons in V1 (24.7%) and 57 extrastriate cells (37.8%). Therefore, we find that roughly 20–25% of the neurons are nonmonotonic in their contrast response overall. Figure 3 shows the monotonicity indices for all nonmonotonic neurons. There is a smooth distribution, with very few neurons being extremely nonmonotonic—only 8 V1 neurons (3.0%) and 7 V2 neurons (3.4%) had a response at full contrast below 50% of their maximum response. 

Figure 3

 

Nonmonotonicity index for neurons in V1 and V2. Nonmonotonicity indices of neurons (with values that are less than 1) are shown. A value of 1 represents a monotonic neuron; a value of 0 represents a neuron whose response at maximum contrast was the same as its response at zero contrast.

Experiment 2: The importance of contrast nonlinearities

Although the saturation of neurons in LGN might be explained by the need to optimize responses to the input contrasts and the need to match that to natural scene statistics (Tadmor & Tolhurst, 2000), this does not account for the shape of the functions found in later visual areas (Clatworthy et al., 2003). In particular, the above finding that 20–25% of primate V1 neurons are nonmonotonic in their contrast response functions would suggest that contrast gain controls are only part of the story for nonlinear contrast responses. This raises the question of what other purpose (or purposes) they might serve. 

An additional potential use for the saturating functions is in detecting the conjunction of certain features. Peirce and Taylor (2006) recently described the existence of mechanisms responding selectively to the presence of a plaid, the conjunction of a pair of overlapping gratings. The mechanism we envisaged for such a detector is similar to that proposed by Olzak and Thomas (1999) and involves a bank of linear, or quasilinear, filters followed by some output nonlinearity (Olzak and Thomas break the saturating nonlinearity into two separate nonlinearities, although for our more limited purposes, this is unnecessary). The nonlinear outputs of the first-level units are then summed by some conjunction detector. 

In keeping with the literature, Olzak and Thomas (1999) suggest that the nonlinearities are important for gain control and normalization. They are, however, also essential for another aspect of the mechanism they describe, namely, for the summing circuit to discriminate between the case of one input channel being stimulated at maximum contrast and that of two input channels each stimulated at 50% (see Figure 4). For this method of discriminating a plaid from one of its components, some form of compressive nonlinearity is essential rather than merely optimal. Furthermore, the neurons in the input layer of the model (e.g., V1 cells) are likely to reduce their response when a second stimulus is present, through the effect of cross-orientation inhibition (e.g., Bonds, 1989; Freeman, Durand, Kiper, & Carandini, 2002; Morrone, Burr, & Maffei, 1982), which exacerbates the problem. These factors can be overcome by appropriate compressive nonlinearities on the outputs of the neurons potentially like those typically found in V1. 

Figure 4

 

The effects of saturating (compressive) nonlinearities on the inputs to a conjunction detector. (A) A linear summing mechanism stimulated at 50% in each input channel will respond equally to stimulation at 100% in either channel, not responding selectively to the conjunction. (B) Adding a nonlinearity to the input channels, as found in V1 neurons, causes the conjunction to provide greater stimulation than either component is capable of alone.

Methods

In considering the extent to which the nonlinearities observed in real neurons might be useful in the selective detection of conjunctions, a further index (conjunction selectivity index [CSI]) was generated. This measures how much more strongly a summing circuit would respond to a full-contrast plaid (for which the neuron being tested provided one of the input channels) versus a single full-contrast grating (for which this neuron is the entire input). This index is naturally defined as  

$\mathrm{C}\mathrm{S}\mathrm{I}=\frac{2R_{50}-R_{100}}{R_{100}},$

(7)

where R ₅₀ represents the response to a grating of 50% contrast and R ₁₀₀ represents the response to a full-contrast grating. In both cases, the response is taken as either the amplitude of modulation at the temporal frequency of stimulation or the spike rate over and above the baseline rate of the neuron, whichever is greater. The index is bounded by zero and positive infinity. A CSI of 0 represents a linear neuron, for which a subsequent summing circuit would respond identically to a plaid and grating; a CSI of 1 indicates perfect saturation in the contrast range 50–100% such that the summing circuit's response to the plaid is 100% greater than the response to a grating. A CSI greater than 1 indicates a supersaturating cell, capable of providing inputs to a very selective summing circuit. 

It should be noted that the CSI, as described, makes the assumption that the presence of the second grating in the stimulus does not affect the response of the neuron, which we know to be untrue in many neurons. In particular, cross-orientation inhibition causes many cells to reduce their response to the preferred grating in the presence of an additional grating stimulus (e.g., Bonds, 1989; Freeman et al., 2002; Morrone et al., 1982). The CSI values here represent the upper limit of the potential conjunction selectivity of the system in such cells. For some cells at least, the effect of the saturation is enough, however, to counter this effect of cross-orientation suppression. Indeed, for the sample cell shown by Freeman et al. (2002, their Figure 2A), the responses to the plaid and the full-contrast grating were 42 and 65 ips, respectively, giving a “true” CSI of 0.30, even in the presence of suppression. 

Results

The extent to which neurons might be useful inputs to conjunction-selective mechanisms is remarkably variable. Figure 5 shows the CSI for V1 simple cells (defined as those whose amplitude of modulation is greater than their maintained spike rate during stimulation), V1 complex cells, and V2 cells. The index represents the extent to which a summing circuit, using this neuron as an input, would be able to discriminate between two weakly stimulated channels and one strongly stimulated channel. The three classes of neuron seem to differ very little in this respect. A number of neurons (80 of the 541 cells studied) show a negative CSI, indicating that rather than saturating, the neuron is actually showing an expansive nonlinearity in the contrast range 0.5–1.0. The majority, however, were nonlinear in a manner useful for conjunction-selective mechanisms. The median CSI was 0.50, such that a full-contrast, two-component plaid might provide a subsequent circuit with 50% more stimulation than a full-contrast grating. Ninety neurons (16%) supersaturated, such that they would provide more than 100% extra summed output during presentation of a plaid. 

Figure 5

 

The CSI for V1 simple-cell neurons (top), V1 complex cells (middle), and V2 cells (bottom). Note that the populations differed very little in this index, all showing a median CSI of roughly 0.5.

General discussion

Previous studies have shown that supersaturation, although not often discussed, is a surprisingly common occurrence in the cat visual cortex (Ledgeway et al., 2005; Li & Creutzfeldt, 1984). The current study shows that it is also very common in macaque visual cortex. There are several reasons why the nonmonotonicity has been largely ignored. The first is simply that extremely nonmonotonic neurons are rare, and minor deviations from monotonicity might occasionally result from noise in the measurements. The second is that high contrasts are sampled sparsely, if at all, by most laboratories. The third is that supersaturation is not captured by the standard function used to represent contrast response functions. The current study adds a second exponent parameter to the traditional Naka–Rushton model, which improves model fits for nearly all neurons and is more efficient, in terms of fit quality normalized by degrees of freedom, for roughly half the neurons studied. Critically, it captures the fact that neurons do not necessarily show monotonically increasing responses to contrast. Although the traditional model serves as a reasonable first approximation to the responses for many cells, its inability to supersaturate leads necessarily to failures in characterizing the 20–25% of neurons found to show this behavior. 

A final reason is that we have had no real explanation for the existence of such a response. Under the dominant hypothesis that the contrast response function serves to optimize the dynamic response range of the neuron, the existence of cells that supersaturate is something of a mystery and might reflect an inappropriately strong contribution of the normalization pool. If this were the case, then it seems surprising that such a large percentage of the neurons supersaturate. If, however, we consider the saturating characteristic of the contrast response as having an additional function in generating selectivity for the conjunction of certain features, then the existence of neurons that supersaturate is less problematic; such cells would be more selective in their detection of compounds. 

This coding mechanism is not necessary for all feature combinations, only those for which the inputs are not independent (it might be considered advantageous in other circumstances, but not necessary). For instance, where the features are spatially overlapping gratings (plaids), the contrast of each component is limited by the contrast of the other, and a saturating nonlinearity is critical in selectively detecting the compound. For features that are independent (e.g., occurring in different eyes or different parts of the retina), a simple linear combination might suffice in detecting the compound, at least for high net contrasts. The expansive, saturating, and supersaturating contrast response curves might therefore be viewed as arising from neurons involved in different forms of processing or different types of feature rather than as suboptimal, optimal, and erroneous forms of contrast gain control mechanisms, respectively. This, of course, supposes that conjunction-selective neurons do exist. Although the psychophysical evidence suggests that such mechanisms are present in the human visual system (Peirce & Taylor, 2006), the author is not aware of any electrophysiological evidence as yet for plaid, form-selective neurons (which respond to a particular plaid more than to either of its isolated components). 

This is also not the only way to perform a nonlinear conjunction-selective combination of signals. An obvious alternative is to combine signals multiplicatively rather than additively (Torre & Poggio, 1978). Computationally, this seems a very different proposal, but its biological implementation may actually involve the summation of two saturating (specifically log-transformed) signals (Tal & Schwartz, 1997). Potentially, these transformations could be performed within a neuron on the postsynaptic potentials, and there is some evidence for this in certain neurons in the fly (Gabbiani, Krapp, Koch, & Laurent, 2002), owl (Pena & Konishi, 2001), and rabbit (Taylor, He, Levick, & Vaney, 2000). No evidence has been shown for such transformations in the primate, however, and it may not be necessary because the outputs of neurons in the previous stages of processing have already endowed signals with a compressive nonlinearity. 

Neurons in the LGN show less extreme saturation (especially the P cells) and do not supersaturate. Their contrast response function seems to be well suited to the optimal coding of natural contrasts (Tadmor & Tolhurst, 2000). This difference with cortical neurons presumably reflects some difference in the role of the contrast response function for cortical and subcortical neurons. It may, for example, reflect a greater degree of independence in stimulation for LGN neurons, resulting from less spatial overlap between receptive fields than found for the larger receptive fields of the cortex. 

The utility of the nonlinear combination of certain signals is not a new concept; computational-modeling and image-processing techniques have frequently used these methods. In visual neuroscience, however, their importance has been somewhat obscured by the consideration of nonlinearities for gain-control mechanisms. Although we can never know the true computational purpose of such signal transformations, this study points to the fact that contrast gain control is not likely to be their sole purpose. It also highlights the problem of using, as a model to characterize the contrast response function, an equation that fundamentally fails to account for the response characteristics of a large proportion of the cortical visual neurons we encounter. 

Acknowledgments

This research was supported by a grant from the BBSRC, UK (BB/C50289X/1), United Kingdom. Data were collected in Peter Lennie's laboratory at NYU, with grants from the NIH, USA (EY04440 and EY13079). Many thanks to J. Forte, J. Kraft, J. Krauskopf, S. Solomon, C. Tailby, and B. Webb for allowing the neuronal data to be used in this analysis. Also, thanks to B. Webb and T. Ledgeway for their comments. 

Commercial relationships: none. 

Corresponding author: Jonathan W. Peirce. 

Email: [e-mail]. 

Address: Nottingham Visual Neuroscience, School of Psychology, University of Nottingham, Nottingham, NG7 2RD, United Kingdom. 

Footnote

Footnotes

¹  In fact, the equation given by Naka and Rushton did not include the exponent n and was therefore identical to that of Michaelis and Menton. Equation 1 is, however, used by most physiology laboratories and still typically referred to as the Naka–Rushton equation.