## Abstract

Connectomes abound, but few for the human spinal cord. Using anatomical data in the literature, we constructed a draft connectivity map of the human spinal cord connectome, providing a template for the many calibrations of specialized behavior to be overlaid on it and the basis for an initial computational model.

A thorough literature review gleaned cell types, connectivity, and connection strength indications. Where human data were not available, we selected species that have been studied. Cadaveric spinal cord measurements, cross-sectional histology images, and cytoarchitectural data regarding cell size and density served as the starting point for estimating numbers of neurons. Simulations were run using neural circuitry simulation software.

The model contains the neural circuitry in all ten Rexed laminae with intralaminar, interlaminar, and intersegmental connections, as well as ascending and descending brain connections and estimated neuron counts for various cell types in every lamina of all 31 segments. We noted the presence of highly interconnected complex networks exhibiting several orders of recurrence. The model was used to perform a detailed study of spinal cord stimulation for analgesia.

This model is a starting point for workers to develop and test hypotheses across an array of biomedical applications focused on the spinal cord. Each such model requires additional calibrations to constrain its output to verifiable predictions. Future work will include simulating additional segments and expanding the research uses of the model.

## 1  Introduction

Connectomes increasingly abound, but few yet for the human spinal cord and none of its large-scale structure, although recent work (e.g., to crystallize neuron types with modern methods and create a genomic-based neuronal atlas in the mouse; Haring et al., 2018), has begun this effort despite its critical and complex role in processing sensory input, motor control, and repertoire of ancient neural motifs such as clocking and pattern generation (Modha & Singh, 2010; Borisyuk, Al Azad, Conte, Roberts, & Soffe, 2011; Marc, Jones, Lauritzen, Watt, & Anderson, 2012; Stobb, Peterson, Mazzag, & Gahtan, 2012; Burns et al., 2013; Kandel, 2013; Blumenfeld, Bliss, Perez, & D'Esposito, 2014; Roberts et al., 2014; Kuan et al., 2015; Meinertzhagen, 2016; Ryan, Lu, & Meinertzhagen, 2016; Barch, 2017; Eichler et al., 2017; Lerman-Sinkoff et al., 2017; Noori et al., 2017). The overall cytoarchitecture and fine-scale physiology have been described (Nógrádi & Vrbová, 2013), and we are learning increasingly more about the circuits, cell membrane properties, numbers and dimensions of axonal branches and collaterals, and synaptic locations (e.g., proximal versus distal) specific to the different classes of spinal cord neurons (Remahl, Cullheim, & Ulfhake, 1985; Enjin et al., 2012; Zhang et al., 2015; Smith et al., 2016; Koch, Acton, & Goulding, 2017). Computer modeling is helping to illuminate likely circuit architecture (Danner, Shevtsova, Frigon, & Rybak, 2017). These data permit a sufficient level of detail from which to construct an initial large-scale human connectome model, similar to ones made for other species, to act as a starting point template for local small-scale connectivity maps and increasingly refined large-scale models.

An early perspective suggested that connectomes fall into two categories of granularity: macro-connectomes, based on imaging data at the highest resolution permitted by technology (e.g. 1–3 mm voxels from resting state and diffusion functional MRI, which connect gray matter regions on that scale) and micro-connectomes at scales of less than 1 mm$3$ (Van Essen et al., 2013). Yet there are hybrid approaches, such as the Allen Brain Project (Kuan et al., 2015) and an increasing number of specialized models at different scales appearing in the literature into which ours falls, driven by the goal of pragmatically modeling neurological disorders and their treatment with available data (Siekmeier & Vanmaanen, 2013; Kucyi & Davis, 2016; Earley, Uhl, Clemens, & Ferre, 2017; Sharma et al., 2017; Tinaz, Lauro, Ghosh, Lungu, & Horovitz, 2017; Van Der Horn et al., 2017). The most audacious connectome vision is whole brain emulation (Markram, 2006; Sandberg & Bostrom, 2008).

Modeling challenges include the typical one of deciding the level of detail that captures the phenomena relevant to testing a particular hypothesis or answering a given research question and is practical computationally (Arle, Shils, & Mei, 2008).

Once the spinal cord connectome was compiled from the literature and various calculations made (see section 2), we implemented it in our Universal Neural Circuitry Simulator (UNCuS; Arle, 1992; Arle et al., 2008; Arle, Carlson, Mei, Iftimia, & Shils, 2014). This software is designed to provide a reasonable compromise between biological realism and computational efficiency for clinical applications. Using this software, simulations researchers can perform in a practical time frame on a standard laptop computer. In UNCuS, neurons are based on a modified MacGregor spiking neuron model (Macgregor, 1987, 1993), which simplifies the classic Hodgkin-Huxley model for efficient computation at the neural circuit level (Arle et al., 2008). UNCuS has been used, for example, to model the circuitry of the cochlear nucleus (Arle & Kim, 1991; Arle, 1992), the basal ganglia in Parkinson's patients undergoing deep brain stimulation treatment (Arle & Shils, 2008), pain signaling in the dorsal horn and its mitigation via spinal cord stimulation (Arle et al., 2014), and the effects of drugs on depressive disorder (Arle & Carlson, 2016).

As an aid in creating the spinal cord connectome and helping to acquire an essential overview of the connectivity, we created visual schematics of the circuitry in PowerPoint slides (see presentation 3 in the online supplements). We employed a bottom-up rather than top-down modeling approach, only incorporating the established connectivity and low-level salient features of the neural circuitry, so that the model would be agnostic to any given theory of motor or sensory function, such as the gate control theory of pain (Melzack & Wall, 1965; Arle et al., 2014; Linderoth & Foreman, 2017), in an attempt to let higher-order phenomena emerge of their own accord. The model presented here contains about 500 distinct neural and fiber groups and 11,000 connections, on the order of the IBM CoCoMac connectome of the macaque brain, which contains 383 neural groups and 6602 connections (Modha & Singh, 2010), scaled 1:50 to reduce computational resources required to run the model, in UNCuS yields 335,000 neurons.

## 2  Methods

We present here a synopsis of our methodology for constructing the spinal cord model (greater detail can be found in the online supplement and Iftimia, 2013). Constructing the connectome consists of these steps, with examples given in our tables and full data sets in the online supplement:

1. Estimating the volume of each spinal cord segment (presentations 1 and 2)

2. Estimating the percent by volume of gray matter in each segment (presentation 2)

3. Determining the segment distribution of neurons (presentation 2)

4. Determining the distribution of neurons in the Rexed lamina in each segment (presentation 2)

5. Determining cell types and connectivity throughout the spinal cord (datasheet 1)

6. Constructing a neural circuitry model of the spinal cord (datasheets 1, 2, 3, 4; presentation 3)

### 2.1  Estimating the Volume of Each Spinal Cord Segment

Cadaveric measurements of spinal cord segment dimensions can be found in the literature (Ko, Park, Shin, & Baek, 2004). This study used a relatively small sample size of 15 cadavers, all of Korean in origin, the only study of individual cord segment dimensions that we found. In the Ko study, cross-sectional cuts were made along the spinal cord, and transverse and sagittal measurements taken at the imaginary lines delineating the proximal (rostral) end of every individual segment, using the spinal nerves associated with each segment as reference points. The rostro-caudal lengths of the segments were also tabulated.

Segments C1 and C2 (first and second cervical) and CC (coccygeal) were not measured in the Ko study, so we extrapolated their dimensions in Microsoft Excel using a best-fit polynomial curve for each independent axis. For segment CC, we also considered that it reduces to 2 mm at its distal (caudal) apex (Drake & Gray, 2008). Clinical experience and studies of CSF thickness (Holsheimer, Barolat, Struijk, & He, 1995; Holsheimer & Barolat, 1998) lead us to note as a limitation that if Koreans are, on average, smaller than a given other population group, the dimensions of their spinal cord may be smaller as well, but we know of no evidence that a smaller individual has less white or gray matter than a larger person, although this is the case with certain disorders such as multiple sclerosis (Liu, Edwards, Gong, Roberts, & Blumhardt, 1999).

Presentation 1 describes the mathematical basis of these results. Briefly, to improve accuracy, we approximated each cross-section as an ellipse and each segment as a truncated cone rather than as an elliptical cylinder. Presentation 2 gives the set of volume estimates for all segments of the cord and Table 1 shows an example for the cervical region. (Note that all presentations are in the online supplement. See the appendix for a list of these presentations.)

Table 1:
Spinal Cord Segment Volumes (mm$3$) from Cross-Sectional and Length Estimate (mm).
 Number of Segment Transverse Sagittal Length Area Volume Neurons C1 11.3 7.7 14.3 68.3 994.3 5,833,495 C2 11.4 7.9 14.8 70.7 1053.8 7,098,515 C3 11.7 7.8 15.3 71.7 1112.7 7,978,838 C4 12.2 7.7 14.0 73.8 1034.2 7,415,938 C5 12.9 7.3 15.5 74.0 1112.1 7,732,883 C6 12.3 7.2 13.2 69.6 884.5 5,573,698 C7 11.9 6.9 12.7 64.5 789.2 4,458,698 C8 11.2 6.8 13.7 59.8 740.5 3,700,842
 Number of Segment Transverse Sagittal Length Area Volume Neurons C1 11.3 7.7 14.3 68.3 994.3 5,833,495 C2 11.4 7.9 14.8 70.7 1053.8 7,098,515 C3 11.7 7.8 15.3 71.7 1112.7 7,978,838 C4 12.2 7.7 14.0 73.8 1034.2 7,415,938 C5 12.9 7.3 15.5 74.0 1112.1 7,732,883 C6 12.3 7.2 13.2 69.6 884.5 5,573,698 C7 11.9 6.9 12.7 64.5 789.2 4,458,698 C8 11.2 6.8 13.7 59.8 740.5 3,700,842

Note: Cervical region estimates for the entire spinal cord are given in presentation 2.

It might be thought that some method of modern imaging could produce more accurate data; to our knowledge, such has not been done and is less feasible than may be imagined (Wheeler-Kingshott et al., 2014).

### 2.2  Estimating the Volume of Gray Matter throughout the Spinal Cord

We used two authoritative secondary sources, Netter (2014) and Carpenter (1996) corroborated by several other sources (Greenstein & Greenstein, 2000; Noback, 2005; Drake & Gray, 2008), to ascertain the percentage of gray matter throughout the spinal cord. These sources illustrate cross-sections of the cord at various levels in diagrammatic form based on experimental data or directly as scanned histology slides. We used public domain image-processing software (ImageJ, http://rsb.info.nih.gov), to determine the area of the gray matter and the total area in each image, allowing us to calculate the percentage.

Since we approximated each spinal cord segment as a truncated elliptical cylinder (see section 3.1 and presentation 1), that methodology carried through to the gray matter volume estimates, which resulted from multiplying the percentage gray matter area that we calculated for each cross-sectional level of the cord times the segment volume.

Not all sources depicted all cord levels, so we extrapolated using a best-fit polynomial curve in Excel to fill missing data and averaged values when there was a significant difference reported across sources for the same level. The following trend was observed: the percentage of gray matter is not constant, so the area of gray matter does not simply vary proportionally with cross-sectional area. Rather, the proportion increases in the thick cervical and lumbosacral enlargements and decreases in the thinner thoracic region.

Presentation 2 gives the complete data set for each segment, and an example is shown in Table 2 for the thoracic region.

Table 2:
Gray Matter Volume Estimates: Thoracic Region.
 Segment Percent Gray Gray Matter Segment Volume Matter Volume T1 734.2 20 146.8 T2 679.9 17 115.6 T3 718.7 15 107.8 T4 755.2 14 105.7 T5 726.0 13 94.4 T6 750.0 13 97.5 T7 717.3 14 100.4 T8 647.7 16 103.6 T9 663.9 18 119.5 T10 597.1 21 125.4 T11 525.4 24 126.1 T12 516.1 27 139.3
 Segment Percent Gray Gray Matter Segment Volume Matter Volume T1 734.2 20 146.8 T2 679.9 17 115.6 T3 718.7 15 107.8 T4 755.2 14 105.7 T5 726.0 13 94.4 T6 750.0 13 97.5 T7 717.3 14 100.4 T8 647.7 16 103.6 T9 663.9 18 119.5 T10 597.1 21 125.4 T11 525.4 24 126.1 T12 516.1 27 139.3

### 2.3  Determining the Segment Distribution of Neurons in the Spinal Cord

It is estimated that there are on the order of 10$8$ neurons in the spinal cord (Longstaff, 2000) constituting its gray matter, while its white matter is composed of axons. If we assume that the average gray matter density from one segment to another is constant within each Rexed lamina, then the number of neurons in any segment depends solely on the volume of its gray matter. Multiplying the percentage of or per segment by the total number of neurons gives the counts by segment in presentation 2 (with an example for the thoracic region given in Table 5).

### 2.4  Determining the Laminar Distribution of Neurons in the Spinal Cord

There are significant differences in neuron size and density among laminae at the same cord level, which are relatively constant along its length. To take these variations into consideration, we constructed a set of quantitative guidelines by which to weight the laminar neuron counts (see Tables 3 and 4; Brown, 1981; Nógrádi & Vrbová, 2013). Multiplying the neuron size factor by the neuron density factor from Table 3 gives the overall weighting factors in Table 4. These weighting factors should be regarded as semiquantitative and to be refined in the future.

Table 3:
Factors Used in Estimating Numbers of Neurons in Each Rexed Lamina.
 Neuron Size Neuron Density Large (50–100 $μ$m) 1/10 High 2 Medium-large (20–50 $μ$m) 1/3 Medium-high 3/2 Medium (10–20 $μ$m) 1 Medium 1 Medium-small (5–10 $μ$m) 3 Medium-low 2/3 Small (1–5 $μ$m) 10 Low 1/2
 Neuron Size Neuron Density Large (50–100 $μ$m) 1/10 High 2 Medium-large (20–50 $μ$m) 1/3 Medium-high 3/2 Medium (10–20 $μ$m) 1 Medium 1 Medium-small (5–10 $μ$m) 3 Medium-low 2/3 Small (1–5 $μ$m) 10 Low 1/2
Table 4:
Weighting Guidelines for Estimating Laminar Neuron Counts.
 Average Average Weighting Lamina Neuron Size Neuron Density Factor I Large Low 1/20 II Small High 20 III Medium Medium-high 3/2 IV Medium-large Medium 1/3 V Medium Medium 1 VI Medium-small Medium 3 VII Medium Medium 1 VIII Medium-small Medium 3 IX Medium-large Medium 1/3 X Medium-small High 6
 Average Average Weighting Lamina Neuron Size Neuron Density Factor I Large Low 1/20 II Small High 20 III Medium Medium-high 3/2 IV Medium-large Medium 1/3 V Medium Medium 1 VI Medium-small Medium 3 VII Medium Medium 1 VIII Medium-small Medium 3 IX Medium-large Medium 1/3 X Medium-small High 6

Multiplying the laminar neuron counts by the appropriate weighting factors changes the total number of neurons, so we normalized the values to maintain the correct total for each segment. All calculations were performed using Excel, allowing for easy adjustments as more accurate measurements appear in the literature. Neuron counts for every lamina of all 31 segments are given in presentation 2, with an example for segment T9 in the thoracic region in Table 5.

Table 5:
Neuron Counts for Segment T9.
 Number of Lamina % Volume Neurons Weighted Scaled I 3 77,901 1,167 23 II 10 259,671 1,556,545 31,131 III 9 233,704 105,067 2,101 IV 10 259,671 25,942 519 V 4 103,868 31,131 623 VI 0 — — — VII 34 882,880 264,613 5,292 VIII 13 337,572 303,526 6,071 IX 11 285,638 28,537 571 X 6 155,802 280,178 5604 Total 100 2,596,706 2,596,706 51,934 Scaling factor 0.020
 Number of Lamina % Volume Neurons Weighted Scaled I 3 77,901 1,167 23 II 10 259,671 1,556,545 31,131 III 9 233,704 105,067 2,101 IV 10 259,671 25,942 519 V 4 103,868 31,131 623 VI 0 — — — VII 34 882,880 264,613 5,292 VIII 13 337,572 303,526 6,071 IX 11 285,638 28,537 571 X 6 155,802 280,178 5604 Total 100 2,596,706 2,596,706 51,934 Scaling factor 0.020

Note: There is no lamina VI in T9.

As is common in modeling, there is a trade-off between actual numbers of neurons in the spinal cord and the computational resources available for simulation, in the context of any given goal. Accordingly, for our initial purposes, we found scaling the model down by a factor of 50 to be desirable (see Table 5, last column).

### 2.5  Determining Cell Types and Connectivity throughout the Spinal Cord

The laminar neuron count results were subsequently scaled down by a factor of 50 for implementation in UNCuS (see Table 5) to conserve computational resources and maintain a reasonable simulation speed. The stellate neuron, which features an approximately linear V-I curve through most of its domain and compressed at domain minimum and maximum (a typical S curve), and whose firing behavior can be modified more easily than nonlinear cell types, was used as a generic default cell type with inhibitory or excitatory polarity. Our emphasis was on calculating cell numbers within groups rather than determining cell types with differing properties, and thus this is a limitation. Classification of cell types in the gray matter has evolved over time and, only very recently perhaps have methods to create an immutable classification emerged (Willis & Coggeshall, 2004b; Usoskin et al., 2015).

Numbers and types of afferent fibers, with an emphasis on those involved in spinal cord stimulation for treatment of neuropathic pain, were also calculated and are presented in Tables 6 to 8. The spinal cord processes afferent fiber sensory information, some locally in the cord and some in higher brain centers, and sends motor control information via efferent fibers. In modeling the spinal cord and in modeling neuromodulation or “electroceuticals” (modification of the nervous system via electromagnetic fields) generally, one must acquire the best information available of fiber function (i.e., type of sensation or motor signal carried), diameter (from which conduction velocity and activation and blocking thresholds may be estimated), location within a region or fascicle (from which to estimate recruitment via finite element modeling of the electric field through tissue and active fiber models to compute threshold), numbers (to estimate total driving force), and their sources (peripheral receptors or interneurons) and targets, to uncover efficacious and inefficacious (side effects) signaling (see Presentation 8 and Carlson et al. (2018)).

Table 6:
Numbers and Types of Large Diameter Afferent Fibers in the Spinal Cord.
 Fiber Type/Sensation Number of Collaterals Conduction Speed/Diameter Rexed Laminae Innervated/Bouton Density Ia afferents (1) Stretch, length, or rate of length change of muscle 3–9, mean 5.8 $±$1.9 $>$80 m/s $>$13 $μ$m VI, VII, IX Ib afferents (2) Contraction or stretch of muscle 5–11, mean 7.4 $±$2 $>$80 m/s $>$13 $μ$m IV, V, VI, dorsal VII 56–384, strongest in medial VI FAII Pacinian corpuscles (5) Tickling or vibration over 20–50 Hz 4–12, mean 7.11 $±$3.55 57–75 m/s 9.5–12.5 $μ$m 11.5 and 12.8 in our model 2/3 in III (strong), IV (weaker) 1/3 in V, VI 47–367 FAI rapidly adapting (RA, Krause's end bulb, Meissner corpuscle) afferents (6) Tapping at 1 Hz; flutter at 10 Hz; vibration at 50 Hz Mean 5 54–60 m/s 11.7–13.0 $μ$m 11.5 and 12.8 $μ$m in our model Almost all in III, rest in IV (number of boutons not given) SAII afferents (Ruffini; subtype: claw base) (5) Usually no sensation (occasionally position) Mean 9.89 28–51 m/s 6.1–11.1 $μ$m 8.7 $μ$m in our model III–dorsal IV (40–50) V, VI $<$ ventral IV $<$ III, dIV V–VI (10–20) SAI afferents (Merkel cell of glabrous and hairy skin origin) (3) Sustained pressure over 5–10 Hz 4–16, mean 8.6 Not given in (3); from (4): 33–95 m/s 7.2–15.8 $μ$m 7.5 $μ$m in our model III, IV, dorsal V hairy (16–46, mean 24.2) glabrous (9–26, mean 15.3) Hair follicles (7) Hair movement 2–12, mean 4.8 Three types of low-threshold mechanoreceptors (RA, hair follicles, SAI) lumbar (8) See sensations above RA: 10.06 $±$0.93 hair follicle: 9.33 $±$0.56 SAI: 9.75 $±$0.75 A$β$ hair follicle afferents in different lumbar spinal levels corresponding to 5 dermatomes (9) Hair movement 9.83 $±$0.64 (over 5 dermatomes) A$β$ RA afferents from various parts of the rat foot (10) Tapping at 1 Hz; flutter at 10 Hz; vibration at 50 Hz 10.82–14.00 (range of means)
 Fiber Type/Sensation Number of Collaterals Conduction Speed/Diameter Rexed Laminae Innervated/Bouton Density Ia afferents (1) Stretch, length, or rate of length change of muscle 3–9, mean 5.8 $±$1.9 $>$80 m/s $>$13 $μ$m VI, VII, IX Ib afferents (2) Contraction or stretch of muscle 5–11, mean 7.4 $±$2 $>$80 m/s $>$13 $μ$m IV, V, VI, dorsal VII 56–384, strongest in medial VI FAII Pacinian corpuscles (5) Tickling or vibration over 20–50 Hz 4–12, mean 7.11 $±$3.55 57–75 m/s 9.5–12.5 $μ$m 11.5 and 12.8 in our model 2/3 in III (strong), IV (weaker) 1/3 in V, VI 47–367 FAI rapidly adapting (RA, Krause's end bulb, Meissner corpuscle) afferents (6) Tapping at 1 Hz; flutter at 10 Hz; vibration at 50 Hz Mean 5 54–60 m/s 11.7–13.0 $μ$m 11.5 and 12.8 $μ$m in our model Almost all in III, rest in IV (number of boutons not given) SAII afferents (Ruffini; subtype: claw base) (5) Usually no sensation (occasionally position) Mean 9.89 28–51 m/s 6.1–11.1 $μ$m 8.7 $μ$m in our model III–dorsal IV (40–50) V, VI $<$ ventral IV $<$ III, dIV V–VI (10–20) SAI afferents (Merkel cell of glabrous and hairy skin origin) (3) Sustained pressure over 5–10 Hz 4–16, mean 8.6 Not given in (3); from (4): 33–95 m/s 7.2–15.8 $μ$m 7.5 $μ$m in our model III, IV, dorsal V hairy (16–46, mean 24.2) glabrous (9–26, mean 15.3) Hair follicles (7) Hair movement 2–12, mean 4.8 Three types of low-threshold mechanoreceptors (RA, hair follicles, SAI) lumbar (8) See sensations above RA: 10.06 $±$0.93 hair follicle: 9.33 $±$0.56 SAI: 9.75 $±$0.75 A$β$ hair follicle afferents in different lumbar spinal levels corresponding to 5 dermatomes (9) Hair movement 9.83 $±$0.64 (over 5 dermatomes) A$β$ RA afferents from various parts of the rat foot (10) Tapping at 1 Hz; flutter at 10 Hz; vibration at 50 Hz 10.82–14.00 (range of means)

1. Brown and Fyffe (1978). 2. Brown and Fyffe (1979). 3. Brown and Fyffe (1978). Brown, Rose, and Snow (1978). 4. Willis and Coggeshall (2004a) 5. Brown, Fyffe, Rose, and Snow (1981). 6. Brown, Fyffe, and Noble (1980). 7. Brown, Rose, and Snow (1977). 8. Woolf (1987). 9. Shortland, Woolf, and Fitzgerald (1989). 10. Shortland and Woolf (1993). FAI: Fast-adapting receptor type I (aka rapidly-adapting or RA). RAII: Fast-adapting/rapidly adapting receptor type II (aka Pacinian corpuscle). SA: Slowly adapting receptor. A$β$: A-beta fiber (large-diameter sensory fiber).

For example, in spinal cord stimulation, pain suppression via large-diameter (more than 10 $μ$m diameter) fibers targeting Rexed lamina III or IV inhibitory interneurons, which in turn suppress lamina V wide dynamic range cells that signal pain to higher centers, is believed to be the mechanism of efficacy, while activation of vibration-signaling fast-adapting I and II fibers is believed to cause the undesirable side effect of paresthesia (Holsheimer, 2002; Lee, Bradley, Kormyio, & Moeller-Bertram, 2011; Lee, Hershey, Bradley, & Yearwood, 2011; Arle, Mei, Carlson, & Shils, 2016). By contrast, for neuromodulation of a different part of the nervous system, vagus nerve stimulation for epilepsy, the current belief is that 7–8 $μ$m fibers cause hoarseness (prevalent when the device is on) by inadvertently stimulating the recurrent laryngeal nerve, while 5 $μ$m fibers emanating from the aortic baroreceptor and are relayed via the nucleus tractus solitarius to the locus coeruleus, cause efficacy by stimulating release of epinephrine, and at the ceiling of clinical amplitudes, higher-threshold 3 $μ$m pulmonary efferent fibers cause undesired coughing (Arle, Carlson, & Mei, 2016).

In high-frequency and burst spinal cord stimulation, analgesia is produced while no paresthesia is produced (De Ridder, Plazier, Kamerling, Menovsky, & Vanneste, 2013). The fiber diameters and targets within the spinal cord gray matter in Table 6 allow one to theorize that the recruitment of larger fibers causing paresthesia is blocked while analgesia is caused by recruitment of smaller fibers causing efficacy (Arle, Mei et al., 2016).

In Tables 7 and 8, note the exponential increase in fiber numbers as their diameter decreases. These and other data, including the dermatome mapping of the dorsal column, along with finite-element modeling of the electric field, lead to the theory that relatively few large-diameter fibers (> 10.7 $μ$m) residing in the superficial 0.3 mm of the dorsal column cause efficacy and side effects in low-frequency spinal cord stimulation and that the floor of efficacy is the threshold of a few most lateral large fibers (Holsheimer, 2002; Lee, Bradley et al., 2011; Lee, Harshey et al., 2011; Arle et al., 2014).

Table 7:
A$β$ Fiber Parameters and Distribution, Dorsal Column (DC) Part A.
 Relative Number Collateral Bouton Number of DC Fiber Diameter Incidence of Fibers Per mm$2$ per Fiber per Fiber Terminals Scaled 1:25 5.7 0.25 42,254 7042 7.5 0.1 16,901 2817 8.6 15 2,180,282 87,211 8.7 0.03 5,070 845 9.89 45 2,256,592 90,264 11.5 0.0025 423 70.4 6.1 207 533,535 21,341 12.8 0.0008 135 22.5 6.1 207 170,731 6829 14 0.00025 42 7.04 5.8 220 53,915 2157 Total 0.38355 64,825 10,804 5,195,055 207,802 Estimates DC fibers T10-11 169,014 Scaling factor 25
 Relative Number Collateral Bouton Number of DC Fiber Diameter Incidence of Fibers Per mm$2$ per Fiber per Fiber Terminals Scaled 1:25 5.7 0.25 42,254 7042 7.5 0.1 16,901 2817 8.6 15 2,180,282 87,211 8.7 0.03 5,070 845 9.89 45 2,256,592 90,264 11.5 0.0025 423 70.4 6.1 207 533,535 21,341 12.8 0.0008 135 22.5 6.1 207 170,731 6829 14 0.00025 42 7.04 5.8 220 53,915 2157 Total 0.38355 64,825 10,804 5,195,055 207,802 Estimates DC fibers T10-11 169,014 Scaling factor 25
Table 8:
A$β$ Fiber Parameters and Distribution, Dorsal Column, Part B.
 Fiber Fibers Fibers per Terminals Diam per 300 $μ$m$2$ Terminals Terminals Scaled per $μ$m 300 $μ$m$2$ Dermatome per mm$2$ per 300 $μ$m$2$ 1:25 Dermatome 5.7 2115 176 7.5 846 70 363,380 109,123 4365 9094 8.7 254 21 376,099 112,943 4518 9412 11.5 21 1.8 88,923 26,703 1068 2225 12.8 7 0.6 28,455 8545 342 712 14 2 0.2 8986 2698 108 225 Total 3245 270 865,843 260,013 10,401 21,668
 Fiber Fibers Fibers per Terminals Diam per 300 $μ$m$2$ Terminals Terminals Scaled per $μ$m 300 $μ$m$2$ Dermatome per mm$2$ per 300 $μ$m$2$ 1:25 Dermatome 5.7 2115 176 7.5 846 70 363,380 109,123 4365 9094 8.7 254 21 376,099 112,943 4518 9412 11.5 21 1.8 88,923 26,703 1068 2225 12.8 7 0.6 28,455 8545 342 712 14 2 0.2 8986 2698 108 225 Total 3245 270 865,843 260,013 10,401 21,668

Note: Dermatomes at T9: 12.

To piece together the known wiring details of the spinal cord connectome, we performed a thorough biomedical literature review, obtaining established cell types, connectivity, and other relevant data, to the extent that these important parameters are known. When data were not available on Homo sapiens, we selected the anatomically and phylogenetically closest species that have been studied to date, which should be noted as a limitation of the model. As with most other connectome efforts, much more detail remains to be uncovered in the future. Presentation 3 contains wiring schematics of the overall model and each Rexed lamina that may be helpful to get a visual overview of the model and spinal cord, but to build a local or global model, datasheet 1, which contains all details of the model in a form that can be translated into a project file for UNCuS or other neural circuitry software. Column H of the TopoConnections tab gives a citation for each connection, and the full reference is given in the ReferenceList tab.

### 2.6  Implementation in Neural Circuitry Simulation Software (UNCuS)

The spinal cord connectome contained too many groups and connections to input into our neural circuitry software manually as we had done with previous models via a graphical user interface (GUI). Thus, using the Excel file as a datasheet to specify neural groups, topological and random connectivity, and electrical stimulation, we wrote a program in Java (Oracle Inc.) to translate the Excel file into the UNCuS project file, which is a space- and numerical code-delimited text file (datasheets 2 and 3). Various functions of Excel, in combination with the datasheet format, are very useful in building and running a model; for instance, one can toggle off entire groups or individual connections to isolate a study to a subset of the model or simply to have it run faster (e.g., turn off the ventral horn to do a dorsal horn study, or vice versa). Another useful feature with so many groups and connections is to turn on Excel's filter at the head of each column and filter for subsets of groups or connections. And bold or color-coding connection strengths can highlight model settings.

For the reasons already given concerning the limited information currently available in the literature, the model is approximate in some senses and underconstrained. Second, the “long-range” connectivity establishes the overall map but cannot constrain the model to replicate all spinal cord behaviors; it is a foundation for future work (Modha & Singh, 2010; Blumenfeld et al., 2014). Each behavior (i.e., an input-output specification) must be calibrated with specialized knowledge concerning additional or revised neuron types, connectivity, and connection strength, which are bottom-up calibrations, and empirical measurements of the behavior, which are top-down calibrations.

For example, as an initial sample calibration, we implemented the classic Hoffman reflex, a monosynaptic reflex arc from the Ia inputs through the alpha motor neuron and out to a muscle, used to test corticospinal tract response when the end of the middle finder is flicked and the thumb moves in response (Yap, 1967) (circuitry shown in Figure 1 and calibration in Figure 2. While a simple circuit, calibration involved iterating many combinations of connection strength between the excitatory Ia afferent, inhibitory Renshaw cells and alpha motor neuron. We envision a more streamlined method using an optimization algorithm in the future.

Figure 1:

Schematic of the monosynaptic or Hoffman reflex (H-reflex) highlighted in yellow on our connectome map of ventral horn central pattern generator. The circuit has been known to be monosynaptic for many years due to its 5 ms response time. However, the refractory time—50$+$ ms—for this simple circuit was found only via empirical studies (Yap, 1967). See presentation 3 for a variety of circuitry schematics, including for all Rexed laminae and common motor reflexes that will be the basis of future calibrations.

Figure 1:

Schematic of the monosynaptic or Hoffman reflex (H-reflex) highlighted in yellow on our connectome map of ventral horn central pattern generator. The circuit has been known to be monosynaptic for many years due to its 5 ms response time. However, the refractory time—50$+$ ms—for this simple circuit was found only via empirical studies (Yap, 1967). See presentation 3 for a variety of circuitry schematics, including for all Rexed laminae and common motor reflexes that will be the basis of future calibrations.

Figure 2:

The monosynaptic Hoffman reflex (H-reflex) calibrated in the neural circuitry decerebrate spine model simulation from (Arle et al., 2014). $x$-axis: time (ms). $y$-axis: transmembrane potential (mV). Three postsynaptic action potentials (APs) at 12, 47, and 82 milliseconds (red arrows) were initiated within a subset of cells of the Ia fiber group, representing a somatotopic region. An alpha motoneuron cell fires in response to the Ia fiber APs at 12 ms and 82 ms but does not reach threshold at 47 ms when another AP is initiated in the Ia fiber, reflecting the circuit's reported refractory time of 50$+$ ms (Yap, 1967). Plots showing balanced afferent fiber Ia, Renshaw cell, and alpha motor neuron reported firing rates are given in presentation 3, slides 15–18.

Figure 2:

The monosynaptic Hoffman reflex (H-reflex) calibrated in the neural circuitry decerebrate spine model simulation from (Arle et al., 2014). $x$-axis: time (ms). $y$-axis: transmembrane potential (mV). Three postsynaptic action potentials (APs) at 12, 47, and 82 milliseconds (red arrows) were initiated within a subset of cells of the Ia fiber group, representing a somatotopic region. An alpha motoneuron cell fires in response to the Ia fiber APs at 12 ms and 82 ms but does not reach threshold at 47 ms when another AP is initiated in the Ia fiber, reflecting the circuit's reported refractory time of 50$+$ ms (Yap, 1967). Plots showing balanced afferent fiber Ia, Renshaw cell, and alpha motor neuron reported firing rates are given in presentation 3, slides 15–18.

Further calibrations are required for the model to capture specific spinal cord phenomena, and many calibrations are possible but also necessary depending on the phenomena that simulations intend to capture, such as the motor reflexes and central pattern generators shown schematically in presentation 3, slides 19 to 21 and 25 to 28. The slide schematics should be used only to get an overview of circuits, while the spreadsheet can be used as a starting point for new investigations. And workers will likely be aware of the calibrations required for their local models.

See section 3 and Arle et al. (2014) for an example of a more detailed calibration we performed to investigate the mechanism of action in neuromodulation of the spinal cord for analgesia. Here is more comment on our methods. While not all connectivity is known, a large number of connections can be established from the literature, which constrain the model significantly. There are insufficient empirical data on connection strength between neural groups. One way to overcome this gap is to calibrate the model to replicate reported phenomena while not explicitly causing the model to instantiate a hypothesis or theory; such, if valid, should emerge from disparate, underlying empirical calibrations. A second objective was to minimally alter unknown parameters, such as connection strength, to achieve a given calibration. Third, as in all biological signaling, there is always a noise level in a simple or complex circuit, and a given calibration must exceed it. We emulated a speculated evolutionary process in titrating up one or more connection strengths to minimally exceed background noise to evolve and effect calibration of a known behavior (i.e. input-output specification). A geometrically scaled set of connection strength parameters was helpful in this effort (presentation 3, SC_InputStrength tab). For the neuropathic pain circuit, this largely involved adjusting topographical connection strength between a tonically active lamina II interneuron group, a purported ectopically active lamina IV group targeting the lamina V WDR pain-signaling group, and a large-fiber inhibitory interneuron target in lamina IV (see the supplementary datasheet 1 and see Arle et al., 2014).

We envision an alternative and, in many cases, more efficient process in which this procedure is achieved via a software program. Such a program would take as inputs the desired calibration—for instance, the firing rates of WDR projection neurons in response to nociceptive and non-nociceptive afferent inputs. It would identify the shortest topographic paths between the source and target groups—in this case, A$β$ fibers, their immediate targets in laminae III–IV of the dorsal horn, and their shortest path to lamina V WDR cells, and test incremental increases in connection strength to effect the target firing rate of the WDR cells so as to minimally disturb previous calibrations.

Note that the calibration and validation of the neural circuitry software were a separate, earlier effort than the work presented here (Arle & Kim, 1991; Arle, 1992; Arle et al., 2008). A similar calibration can be performed for many cell types present in the spinal cord or, if a current-voltage relationship is known, it can be input as a lookup table. Using the Excel specification as a starting point, the model may be implemented in other software as well (see datasheets 1–3).

The software outputs up to 15 data files, including electrode stimulation by time step, flags of normal versus hyperactivity, and emitted current by neuron by time step, but the most commonly used output is neural membrane potential time series from which action potential plots can be made (e.g., see Figure 2 and slides in presentation 3). Over 20 parameters have a (pseudo) random component, which in most cases can be toggled on or off by the user, and all may be fixed to ensure consistent repeated simulations by a user-input random number seed. There is also a tab in the datasheet to set up random, as opposed to topgraphical, connections between neural groups.

## 3  Results

The results of our connectome methodology are in the previous tables, with a complete dataset of large-scale spinal cord connectivity in the spreadsheet in supplements. Tables 1 to 5 detail the volume estimates, neuron type weightings, and calculations of numbers of neurons in the Rexed lamina for all spinal cord segments as described in section 2 and the online supplements.

Tables 6 to 8 detail key fiber parameters for the larger (A$β)$ sensory fiber group in the dorsal column, including the type of sensation they carry and estimated number of fibers, including their terminals (Feirabend, Choufoer, Ploeger, Holsheimer, & van Gool, 2002; Lee, Hershey et al., 2011). We gathered information for this group of fibers since they are required for our simulations of neuromodulation for neuropathic pain relief (Arle, Mei et al., 2016). Similar work, for which these tables form an example, would have to be done for other investigations involving the remaining fibers such as those in Willis and Coggeshall (vol. 1, Table 2.2, and text, 2004b).

As an example, we performed a longer and much more involved calibration to investigate how dorsal column stimulation interferes with neuropathic pain signaling (Arle et al., 2014). The calibration began with the circuitry detail initially in the spreadsheet and continued with iterative adjusting of connection weights until the circuit performed some key known aspects and illuminated others about which less was known. We discovered, as others have reported (Willis & Coggeshall, 2004b), that the gate control theory, while illuminative in its basic premises that large sensory fiber stimulation could result in shutting down pain signaling, was incorrect in some key aspects and suggested revised circuit connectivity targeting lamina V WDR cells via large fiber stimulation of laminae III or IV inhibitory interneurons (Arle et al., 2014; Melzack & Wall, 1965; Willis & Coggeshall, 2004b).

In our and similar investigations by others (Holsheimer, 2002; Lee, Bradley et al., 2011; Lee, Harshey et al., 2011), and in an investigation we performed of vagus nerve stimulation for efficacy intended as a model for electroceutical simulations (Arle, Carlson et al., 2016; Carlson et al., 2018), the number of fibers driving the system, their coverage of the dermatomes or other topographically mapped targets, the sensation they provide that may give clues to side effects or their targets, and their diameters that indicate conduction velocity and activation and blocking thresholds are all critical ingredients and constraints to the models. Considerable effort has gone into determining these data by various labs. Some of this information for the spinal cord model, notably that is relevant to our investigations of pain and analgesia circuitry, is provided in Tables 6, 7, and 8.

To get a visual overview of each Rexed lamina's internal and external wiring and an idea of the circuit's complexity, refer to the schematics in presentation 3 (Spinal Cord Schematics—Ventral Horn Calibration).

Table 9 shows a subset of the 474 neural groups, and Table 10 shows a subset of the nearly 11,000 connections between neural groups in the spinal cord connectome that were gleaned from the literature. These are partial glimpses of the full datasheet fields. Readers are encouraged to view the complete sheet in datasheet 1, where a key to abbreviations is provided, as well as the references for the UNCuS neural circuitry software.

Table 9:
Subset of the 474 Neural Groups Contained in the Spinal Cord Connectome Datasheet.
 Group Population Group Neuron Number Start End ID On/Off Group Name Number Number Type of Cells CellID CellID 101 on L1_I_left_I-PYR 1 1 Stellate 80 1 80 102 on L1_I_left_I-WDR 1 2 Stellate 80 81 160 103 on L1_I_left_I-NS 1 3 Stellate 80 161 240 104 on L1_I_left_sI-C-PYR[inh] 1 4 Stellate 80 241 320 105 on L1_I_left_sI-C-WDR[inh] 1 5 Stellate 80 321 400 106 on L1_I_left_sI-C-NS[inh] 1 6 Stellate 80 401 480 107 on L1_I_left_sI-C-PYR[exc] 1 7 Stellate 80 481 560 108 on L1_I_left_sI-C-WDR[exc] 1 8 Stellate 80 561 640 109 on L1_I_left_sI-C-NS[exc] 1 9 Stellate 80 641 720 110 on L1_I_left_C1-PYR 1 10 Stellate 80 721 800 111 on L1_I_left_C2-WDR 1 11 Stellate 80 801 880 112 on L1_I_left_C3-NS 1 12 Stellate 84 881 964 113 on L1_II_left_I-PYR 1 13 Stellate 2960 965 3924 114 on L1_II_left_I-WDR 1 14 Stellate 2960 3925 6884 115 on L1_II_left_I-NS 1 15 Stellate 2960 6885 9844 116 on L1_II_left_sI-C-PYR[inh] 1 16 Stellate 2960 9845 12,804 117 on L1_II_left_sI-C-WDR[inh] 1 17 Stellate 2960 12,805 15,764 118 on L1_II_left_sI-C-NS[inh] 1 18 Stellate 2960 15,765 18,724 119 on L1_II_left_sI-C-PYR[exc] 1 19 Stellate 2960 18,725 21,684 120 on L1_II_left_sI-C-WDR[exc] 1 20 Stellate 2960 21,685 24,644
 Group Population Group Neuron Number Start End ID On/Off Group Name Number Number Type of Cells CellID CellID 101 on L1_I_left_I-PYR 1 1 Stellate 80 1 80 102 on L1_I_left_I-WDR 1 2 Stellate 80 81 160 103 on L1_I_left_I-NS 1 3 Stellate 80 161 240 104 on L1_I_left_sI-C-PYR[inh] 1 4 Stellate 80 241 320 105 on L1_I_left_sI-C-WDR[inh] 1 5 Stellate 80 321 400 106 on L1_I_left_sI-C-NS[inh] 1 6 Stellate 80 401 480 107 on L1_I_left_sI-C-PYR[exc] 1 7 Stellate 80 481 560 108 on L1_I_left_sI-C-WDR[exc] 1 8 Stellate 80 561 640 109 on L1_I_left_sI-C-NS[exc] 1 9 Stellate 80 641 720 110 on L1_I_left_C1-PYR 1 10 Stellate 80 721 800 111 on L1_I_left_C2-WDR 1 11 Stellate 80 801 880 112 on L1_I_left_C3-NS 1 12 Stellate 84 881 964 113 on L1_II_left_I-PYR 1 13 Stellate 2960 965 3924 114 on L1_II_left_I-WDR 1 14 Stellate 2960 3925 6884 115 on L1_II_left_I-NS 1 15 Stellate 2960 6885 9844 116 on L1_II_left_sI-C-PYR[inh] 1 16 Stellate 2960 9845 12,804 117 on L1_II_left_sI-C-WDR[inh] 1 17 Stellate 2960 12,805 15,764 118 on L1_II_left_sI-C-NS[inh] 1 18 Stellate 2960 15,765 18,724 119 on L1_II_left_sI-C-PYR[exc] 1 19 Stellate 2960 18,725 21,684 120 on L1_II_left_sI-C-WDR[exc] 1 20 Stellate 2960 21,685 24,644
Table 10:
Subset of the Nearly 11,000 Topographical Connections between Neural Groups Contained in the Spinal Cord Connectome Datasheet.
 Source Target Connection Strength Axon Dendritic Compartment Source Group Name Target Group Name Minimum Maximum Range Delay C1 C2 C3 C4 C5 C6 C7 C8 C9 C10 L1_I_left_I-LTM L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-LTM L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-WDR L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-WDR L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-NS L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-NS L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-WDR[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-WDR[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-NS[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-NS[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0
 Source Target Connection Strength Axon Dendritic Compartment Source Group Name Target Group Name Minimum Maximum Range Delay C1 C2 C3 C4 C5 C6 C7 C8 C9 C10 L1_I_left_I-LTM L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-LTM L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-WDR L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-WDR L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-NS L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_I-NS L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-WDR[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-WDR[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-NS[inh] L1_II_left_V-SL-WDR[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-NS[inh] L1_II_left_V-SL-NS[exc] 4 6 10 0.25 0 0 0 0 1 0 0 0 0 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-LTM[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-WDR[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-SL-WDR[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-SL-NS[inh] 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-FIL 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0 L1_I_left_C-sI-NS[exc] L1_II_left_V-CUR 4 6 10 0.25 0 0 0 0 0 0 0 0 1 0

Note: The connection strength and dendritic compartment methodologies are explained in the references for the UNCuS neural circuitry software (Arle & Kim, 1991; Arle, 1992; Arle et al., 2008, 2014).

## 4  Discussion

### 4.1  Summary

Dynamic modeling of the spinal cord connectome has innumerable specialized and varied theoretical applications, such as analyzing the mechanisms of cyclic behavior (e.g., limb and gait control in locomotion; Grillner, Wallen, Saitoh, Kozlov, & Robertson, 2008; Danner et al., 2017), and neurological disorders such as neuropathic pain (Holsheimer, 2002; Arle et al., 2014). There are numerous translational aspects, and toward that end we used the model as a starting point to perform further calibration of afferent fiber input, acute peripheral pain, and centrally sensitized neuropathic pain and then explored how dorsal column stimulation ameliorates neuropathic pain to inform the placement and programming of the stimulating electrodes and increase clinical efficacy (Arle et al., 2014).

### 4.2  Limitations

There are several limitations to our model, a work in progress designed to readily incorporate new data that appear in the literature. Notably, the model is intended to be a comprehensive overall map that includes all significant neural groups and connections known at the time of writing but is not represented to be calibrated for all spinal cord behaviors, which is unfeasible for any one lab. That said, such topographical connectivity maps are increasingly proving to be correlated to functionality—the “functional connectome” (Wang, Dai, Gong, Zhou, & He, 2015; Grandjean, Zerbi, Balsters, Wenderoth, & Rudin, 2017; Ferrario, Merrison-Hort, Soffe, & Borisyuk, 2018). A desired spinal cord investigation can take the model, apply calibrations known to specialists to constrain it, and then run sweeps across specific parameters to which dependent variables of interest are believed to be sensitive. In that sense it is relatively easy and straightforward to change field values in the datasheet (see datasheet 1) to effect a calibration. Programming a specific neuron type may be tedious, while using a simplifying assumption, such as a current-voltage lookup table, to effect the same end may be simple.

Here are some general limitations of which we are aware. Although we mapped the ascending and descending loops with the brain (see datasheet 1, BrainConnections tab), we did not include these as active pathways within our computational model. This decision stemmed from the hypothesis that retrograde dorsal column stimulation alone is sufficient to yield pain relief (Linderoth & Meyerson, 2010). We did not find data on connection strength to be plentiful or quantitative, and using what we did find and equal excitatory and inhibitory connection strength for all other connections resulted in a chronically hyperactive condition. In retrospect, it is unsurprising that a model with approximately 8700 excitatory and 2300 inhibitory connections and a great deal of recurrent connectivity should be inherently hyperactive. Accordingly, we adjusted the default connection strengths so that there existed nonhyperactive background activity essentially by decreasing default-level connection strengths of the both excitatory and inhibitory groups until background activity was at an acceptable low level. Against this backdrop, increased connections can be imposed to implement a circuit that effects a higher signal-to-noise ratio than that default. The inherent hyperactivity flagged the likelihood that descending inhibitory connections, stabilized oscillating circuitry (Buzsaki & Draguhn, 2004; Diba, Amarasingham, Mizuseki, & Buzsaki, 2014; Salkoff, Zagha, Yuzgec, & McCormick, 2015; Buzsaki & Christen, 2016) and other constraints are necessary conditions for nonhyperactive, normal activity. Furthermore, innovative techniques for hypothesizing connectivity (Tang, Sun, Toga, Ringman, & Shi, 2017) may be efficacious while workers await more precise data.

In our simulations based on this model, we scaled down neuron counts by a factor of 50, which yielded sufficient detail for our purposes but is nevertheless a compromise in the biological accuracy of the model. We could not scale lower since the smaller neuron groups had too few neurons and they were overdriven by larger groups, so we scaled back up until their output (e.g., firing rates) seemed reasonable. An alternate method would be to adjust connection strength between neural groups, provided the calibration constraints are satisfied. Generally, adjusting connection strength must be done when the relative numbers of neurons are changed, but not necessarily when the ratios are kept fixed. Scaling can be determined by the computational resources and biological accuracy required for the project at hand. While no biochemistry (neurotransmitters and neuromodulators) was included, various biochemical effects can be represented in neural circuitry software, for example, representing inhibitory effects of gamma-aminobutyric acid (GABA) or glycine, or excitatory effects of glutamate, locally via synaptic strength or globally via overall inhibitory or excitatory factors. The default current voltage curve of a neuron was sufficient for capturing its behavior relative to the investigations we performed. In general, we doubt similar modeling efforts will encompass several biological systems levels (Ramirez-Mahaluf, Roxin, Mayberg, & Compte, 2015; Mclean, Shipley, Nickell, Aston-Jones, & Reyher, 1989; Gerstner & Kistler, 2002). For instance, we have developed several active axon models that capture ion channel gate kinetics, but, being computationally intensive themselves, feed the results of these models into separate, higher-level circuitry models as described here (Arle, Mei et al., 2016).

Only one spinal cord segment has been simulated thus far. The connectivity is almost identical for all segments, and we have mapped out the intersegmental connections and calculated the neuron counts for each segment. Expanding simulation to additional segments will require some additional work and computational resources.

### 4.3  The Future of the Spinal Cord Connectome

While we took the bottom-up approach to the limit by incorporating as many data as possible on the number of fibers and neurons, and the connectivity between neural groups, to further impose constraints on such a complex model, more behavior-level calibrations could be implemented. For example, one avenue would combine pain and muscle circuitry in instantiating the flexor-crossed extensor reflex, where painful A-delta and C fiber input triggers lamina V wide dynamic range cells to activate lamina IX motor neurons (e.g., touching a hot stove).

Modeling how motor reflexes are impaired by lesions in various parts of the spinal cord would illuminate the pathophysiology of paralysis and elucidate methods to improve treatment.

Another interesting avenue would include biochemistry to simulate drug action on the spinal cord as a fast and inexpensive precursor to animal and human clinical trials (Arle & Carlson, 2016).

### 4.4  Conclusion: Whole Spinal Cord Emulation

The long-term logical end point of these efforts is whole spinal cord emulation, a sufficiently detailed and complete computational model of the spinal cord connectome that, when run as a dynamic neural simulation, accurately reproduces known sensorimotor functions and other information processing performed by the spinal cord (like Open Worm—Szigeti et al., 2014—or large-scale brain emulation Markram, 2006; Sandberg & Bostrom, 2008). Such a goal can be approached pari passu with whole brain emulation, riding on the same advances in elucidating neurotransmitter function (Ramirez-Jarquin & Tapia, 2018) and micro- and macrocircuit function, for example, via RNA sequencing to classify cell types, two-photon laser scanning and microscopy, high-resolution anatomical magnetic resonance imaging (MRI), T1- and T2-weighted MRI, magneto- and electroencephalography (MEG/EEG), spatial and temporal control via optogenetics, and other new genetic tools, notably to induce gain or loss of function in target circuits; functional, diffusion, resting-state, and task-evoked MRI; probabilistic methods of hypothesizing connectivity; and computer modeling of common circuit motifs such as pattern generators within the coming decades (Deisseroth, 2011; Van Essen et al., 2013; Poldrack & Farah, 2015; Usoskin et al., 2015; Park & Carmel, 2016; Poldrack & Yarkoni, 2016; Danner et al., 2017; Silva, 2017; Tang et al., 2017; Haring et al., 2018; Piatkevich et al., 2018).

As we approach this epochal landmark in computational neuroscience, a vast number of clinical applications and novel branches of research will open up. The model presented here is perhaps the first incipient step toward this ultimate aim.

## Appendix:  Online Supplements

Here are links to a slide deck overview of the project, further detail on our methodology, and supportive information and results. To build a local or global spinal cord model, start from datasheet 1, not the schematic diagrams in presentation 3 or datasheet 4.

•

Presentation 1: Methodology of Spinal Cord Segment Volume Calculations

•

Presentation 2: Spinal Cord Neuron Counts by Segment and Lamina

•

Presentation 3: Slide Deck: Spinal Cord Schematics-Ventral Horn Calibrations

•

Datasheet 1: Complete Spinal Cord Datasheet: Neuron Groups & Topological Connections

•

Datasheet 2: Sample UNCuS Neural Circuitry Datafile

•

Datasheet 3: Format for MS Excel Datasheet of Spinal Cord Model Parameters

•

Datasheet 4: Neuron Key for Spinal Cord Schematics

## Acknowledgments

The work contained in this letter was performed with generous support from the Sydney Family Foundation, Manny Landsman, and the Robert Wise Foundation. Part of the work was performed at the Lahey Clinic, Burlington, Massachusetts. We declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

J.E.A. and J.L.S. conceived of the project and oversaw it. N.I. performed most of the literature research. J.E.A. and L.M. wrote the neural circuitry simulation program. K.W.C. built the model in the neural circuitry software and performed the simulations. K.W.C. and N.I. wrote this letter, which all authors approved.

## References

Arle
,
J. E.
(
1992
).
Neural modeling of the cochlear nucleus.
Ph.D. dissertation, University of Connecticut
.
Arle
,
J. E.
, &
Carlson
,
K. W.
(
2016
).
The use of dynamic computational models of neural circuitry to streamline new drug development
.
Drug Discovery Today: Disease Models
,
19
,
69
75
. doi:https://doi.org/10.1016/j.ddmod.2017.01.002
Arle
,
J. E.
,
Carlson
,
K. W.
, &
Mei
,
L.
(
2016
).
Investigation of mechanisms of vagus nerve stimulation for seizure using finite element modeling
.
Epilepsy Research
,
126
,
109
118
. doi:10.1016/j.eplepsyres.2016.07.009
Arle
,
J. E.
,
Carlson
,
K. W.
,
Mei
,
L.
,
Iftimia
,
N.
, &
Shils
,
J. L.
(
2014
).
Mechanism of dorsal column stimulation to treat neuropathic but not nociceptive pain: Analysis with a computational model
.
Neuromodulation
,
17
(
7
),
642
655
. doi:10.1111/ner.12178
Arle
,
J. E.
, &
Kim
,
D. O.
(
1991
).
Simulations of cochlear nucleus neural circuitry: Excitatory-inhibitory response-area types I–IV
.
Journal of the Acoustical Society of America
,
90
(
6
),
3106
3121
.
Arle
,
J. E.
,
Mei
,
L.
,
Carlson
,
K. W.
, &
Shils
,
J. L.
(
2016
).
High-frequency stimulation of dorsal column axons: Potential underlying mechanism of paresthesia-free neuropathic pain relief
.
Neuromodulation
,
19
(
4
),
385
397
. doi:10.1111/ner.12436
Arle
,
J. E.
, &
Shils
,
J. L.
(
2008
).
Motor cortex stimulation for pain and movement disorders
.
Neurotherapeutics
,
5
(
1
),
37
49
. doi:10.1016/j.nurt.2007.11.004
Arle
,
J. E.
,
Shils
,
J. L.
, &
Mei
,
L. Z.
(
2008
).
Modeling Parkinsonian circuitry and the DBS electrode. I. Biophysical background and software
.
Stereotactic and Functional Neurosurgery
,
86
(
1
),
1
15
.
Barch
,
D. M.
(
2017
).
Resting-state functional connectivity in the human connectome project: Current status and relevance to understanding psychopathology
.
Harvard Review of Psychiatry
,
25
(
5
),
209
217
. doi:10.1097/HRP.0000000000000166
Blumenfeld
,
R. S.
,
Bliss
,
D. P.
,
Perez
,
F.
, &
D'Esposito
,
M.
(
2014
).
Cocotools: Open-source software for building connectomes using the cocomac anatomical database
.
Journal of Cognitive Neuroscience
,
26
(
4
),
722
745
. doi:10.1162/jocn_a_00498
Borisyuk
,
R.
,
Al Azad
,
A. K.
,
Conte
,
D.
,
Roberts
,
A.
, &
Soffe
,
S. R.
(
2011
).
Modeling the connectome of a simple spinal cord
.
Frontiers in Neuroinformatics
,
5
,
20
. doi:10.3389/fninf.2011.00020
Brown
,
A. G.
(
1981
).
Organization in the spinal cord: The anatomy and physiology of identified neurons.
New York
:
Springer-Verlag
.
Brown
,
A. G.
, &
Fyffe
,
R. E.
(
1978
).
The morphology of group Ia afferent fibre collaterals in the spinal cord of the cat
.
Journal of Physiology
,
274
,
111
127
.
Brown
,
A. G.
, &
Fyffe
,
R. E.
(
1979
).
The morphology of group ib afferent fibre collaterals in the spinal cord of the cat
.
Journal of Physiology
,
296
,
215
226
.
Brown
,
A. G.
,
Fyffe
,
R. E.
, &
Noble
,
R.
(
1980
).
Projections from pacinian corpuscles and rapidly adapting mechanoreceptors of glabrous skin to the cat's spinal cord
.
Journal of Physiology
,
307
,
385
400
.
Brown
,
A. G.
,
Fyffe
,
R. E.
,
Rose
,
P. K.
, &
Snow
,
P. J.
(
1981
).
Spinal cord collaterals from axons of type II slowly adapting units in the cat
.
Journal of Physiology
,
316
,
469
480
.
Brown
,
A. G.
,
Rose
,
P. K.
, &
Snow
,
P. J.
(
1977
).
The morphology of hair follicle afferent fibre collaterals in the spinal cord of the cat
.
Journal of Physiology
,
272
(
3
),
779
797
.
Brown
,
A. G.
,
Rose
,
P. K.
, &
Snow
,
P. J.
(
1978
).
Morphology and organization of axon collaterals from afferent fibres of slowly adapting type I units in cat spinal cord
.
Journal of Physiology
,
277
,
15
27
.
Burns
,
R.
,
Roncal
,
W. G.
,
Kleissas
,
D.
,
Lillaney
,
K.
,
Manavalan
,
P.
,
Perlman
,
E.
, …
Vogelstein
,
R. J.
(
2013
).
The Open Connectome Project data cluster: Scalable analysis and vision for high-throughput neuroscience
.
Scientific and Statistical Database Management
. doi:10.1145/2484838.2484870
Buzsaki
,
G.
, &
Christen
,
Y.
(Eds.). (
2016
).
Micro-, meso- and macro-dynamics of the brain.
Cham
:
Springer
.
Buzsaki
,
G.
, &
Draguhn
,
A.
(
2004
).
Neuronal oscillations in cortical networks
.
Science
,
304
(
5679
),
1926
1929
. doi:10.1126/science.1099745
Carlson
,
K. W.
,
Begnaud
,
J.
,
Dokos
,
S.
,
Shils
,
J. L.
,
Mei
,
L.
, &
Arle
,
J. E.
(
2018
).
Simulation as a critical ingredient of the electroceutical paradigm
. In
Proceedings of the 40th Annual International Conference of the IEEE Engineering in Medicine and Biology Society
.
Piscataway, NJ
:
IEEE
.
Carpenter
,
M. B.
(
1996
).
Core text of neuroanatomy.
Baltimore
:
Williams & Wilkins
.
Danner
,
S. M.
,
Shevtsova
,
N. A.
,
Frigon
,
A.
, &
Rybak
,
I. A.
(
2017
).
Computational modeling of spinal circuits controlling limb coordination and gaits in quadrupeds
.
Elife
,
6
. doi:10.7554/eLife.31050
De Ridder
,
D.
,
Plazier
,
M.
,
Kamerling
,
N.
,
Menovsky
,
T.
, &
Vanneste
,
S.
(
2013
).
Burst spinal cord stimulation for limb and back pain
.
World Neurosurgery
,
80
(
5
),
642
649
,
e641
. doi:10.1016/j.wneu.2013.01.040
Deisseroth
,
K.
(
2011
).
Optogenetics
.
Nature Methods
,
8
(
1
),
26
29
. doi:10.1038/nmeth.f.324
Diba
,
K.
,
Amarasingham
,
A.
,
Mizuseki
,
K.
, &
Buzsaki
,
G.
(
2014
).
Millisecond timescale synchrony among hippocampal neurons
.
Journal of Neuroscience
,
34
(
45
),
14984
14994
. doi:10.1523/JNEUROSCI.1091-14.2014
Drake
,
R. L.
, &
Gray
,
H.
(
2008
).
Gray's atlas of anatomy.
Philadelphia
:
Churchill Livingstone Elsevier
.
Earley
,
C. J.
,
Uhl
,
G. R.
,
Clemens
,
S.
, &
Ferre
,
S.
(
2017
).
Connectome and molecular pharmacological differences in the dopaminergic system in restless legs syndrome (Rls): Plastic changes and neuroadaptations that may contribute to augmentation
.
Sleep Medicine
,
31
,
71
77
. doi:10.1016/j.sleep.2016.06.003
Eichler
,
K.
,
Li
,
F.
,
Litwin-Kumar
,
A.
,
Park
,
Y.
,
Andrade
,
I.
,
Schneider-Mizell
,
C. M.
, …
Cardona
,
A.
(
2017
).
The complete connectome of a learning and memory centre in an insect brain
.
Nature
,
548
(
7666
),
175
182
. doi:10.1038/nature23455
Enjin
,
A.
,
Leao
,
K. E.
,
Mikulovic
,
S.
,
Le Merre
,
P.
,
Tourtellotte
,
W. G.
, &
Kullander
,
K.
(
2012
).
Sensorimotor function is modulated by the serotonin receptor 1d, a novel marker for gamma motor neurons
.
Molecular and Cellular Neuroscience
,
49
(
3
),
322
332
. doi:10.1016/j.mcn.2012.01.003
Feirabend
,
H. K. P.
,
Choufoer
,
H.
,
Ploeger
,
S.
,
Holsheimer
,
J.
, &
van Gool
,
J. D.
(
2002
).
Morphometry of human superficial dorsal and dorsolateral column fibres: Significance to spinal cord stimulation
.
Brain
,
125
(
Pt. 5
),
1137
1149
.
Ferrario
,
A.
,
Merrison-Hort
,
R.
,
Soffe
,
S. R.
, &
Borisyuk
,
R.
(
2018
).
Structural and functional properties of a probabilistic model of neuronal connectivity in a simple locomotor network
.
Elife
,
7
. doi:10.7554/eLife.33281
Gerstner
,
W.
, &
Kistler
,
W. M.
(
2002
).
Spiking neuron models: Single neurons, populations, plasticity.
Cambridge
:
Cambridge University Press
.
Grandjean
,
J.
,
Zerbi
,
V.
,
Balsters
,
J. H.
,
Wenderoth
,
N.
, &
Rudin
,
M.
(
2017
).
Structural basis of large-scale functional connectivity in the mouse
.
Journal of Neuroscience
,
37
(
34
),
8092
8101
. doi:10.1523/JNEUROSCI.0438-17.2017
Greenstein
,
B.
, &
Greenstein
,
A.
(
2000
).
Color atlas of neuroscience: Neuroanatomy and neurophysiology.
New York
:
Thieme
.
Grillner
,
S.
,
Wallen
,
P.
,
Saitoh
,
K.
,
Kozlov
,
A.
, &
Robertson
,
B.
(
2008
).
Neural bases of goal-directed locomotion in vertebrates—An overview
.
Brain Research Reviews
,
57
(
1
),
2
12
. doi:10.1016/j.brainresrev.2007.06.027
Haring
,
M.
,
Zeisel
,
A.
,
Hochgerner
,
H.
,
Rinwa
,
P.
,
Jakobsson
,
J. E. T.
,
Lonnerberg
,
P.
, …
Ernfors
,
P.
(
2018
).
Neuronal atlas of the dorsal horn defines its architecture and links sensory input to transcriptional cell types
.
Nature Neuroscience
,
21
(
6
),
869
880
. doi:10.1038/s41593-018-0141-1
Holsheimer
,
J.
(
2002
).
Which neuronal elements are activated directly by spinal cord stimulation?
Neuromodulation
,
5
(
1
),
25
31
.
Holsheimer
,
J.
, &
Barolat
,
G.
(
1998
).
Spinal geometry and paresthesia coverage in spinal cord stimulation
.
Neuromodulation
,
1
(
3
),
129
136
. doi:10.1111/j.1525-1403.1998.tb00006.x
Holsheimer
,
J.
,
Barolat
,
G.
,
Struijk
,
J. J.
, &
He
,
J.
(
1995
).
Significance of the spinal cord position in spinal cord stimulation
.
Acta Neurochirurgica Supplement
,
64
,
119
124
.
Iftimia
,
N. A.
(
2013
).
Spinal cord neural modeling for clinical applications: Senior thesis.
B.A. thesis, Brandeis University
. http://bit.ly/Iftimia_Spinal_cord_neural_modeling
Kandel
,
E. R.
(
2013
).
Principles of neural science.
New York
:
McGraw-Hill Medical
.
Ko
,
H.
,
Park
,
J.
,
Shin
,
Y.
, &
Baek
,
S.
(
2004
).
Gross quantitative measurements of spinal cord segments in human
.
Spinal Cord
,
42
(
1
),
35
40
.
Koch
,
S. C.
,
Acton
,
D.
, &
Goulding
,
M.
(
2017
).
Spinal circuits for touch, pain, and itch
.
Annual Review of Physiology
,
80
. doi:10.1146/annurev-physiol-022516-034303
Kuan
,
L.
,
Li
,
Y.
,
Lau
,
C.
,
Feng
,
D.
,
Bernard
,
A.
,
Sunkin
,
S. M.
, …
Ng
,
L.
(
2015
).
Neuroinformatics of the Allen Mouse Brain Connectivity Atlas
.
Methods
,
73
,
4
17
. doi:10.1016/j.ymeth.2014.12.013
Kucyi
,
A.
, &
Davis
,
K. D.
(
2016
).
The neural code for pain: From single-cell electrophysiology to the dynamic pain connectome
.
Neuroscientist
. doi:10.1177/1073858416667716
Lee
,
D.
,
Bradley
,
K.
,
Kormyio
,
N.
, &
Moeller-Bertram
,
T.
(
2011
).
Investigation of fibers activated in spinal cord stimulation (SCS) using a computational model
. In
Proceedings of the 15th North American Neuromodulation Conference.
Chicago
:
North American Neuromodulation Society
.
Lee
,
D.
,
Hershey
,
B.
,
Bradley
,
K.
, &
Yearwood
,
T.
(
2011
).
Predicted effects of pulse width programming in spinal cord stimulation: A mathematical modeling study
.
Medical and Biological Engineering and Computing
,
49
(
7
),
765
774
. doi:10.1007/s11517-011-0780-9
Lerman-Sinkoff
,
D. B.
,
Sui
,
J.
,
Rachakonda
,
S.
,
Kandala
,
S.
,
Calhoun
,
V. D.
, &
Barch
,
D. M.
(
2017
).
Multimodal neural correlates of cognitive control in the Human Connectome Project
.
NeuroImage
,
163
,
41
54
. doi:10.1016/j.neuroimage.2017.08.081
Linderoth
,
B.
, &
Foreman
,
R. D.
(
2017
).
Conventional and novel spinal stimulation algorithms: Hypothetical mechanisms of action and comments on outcomes
.
Neuromodulation
,
20
(
6
),
525
533
. doi:10.1111/ner.12624
Linderoth
,
B.
, &
Meyerson
,
B. A.
(
2010
).
Spinal cord stimulation: Exploration of the physiological basis of a widely used therapy
.
Anesthesiology
,
113
(
6
),
1265
1267
. doi:10.1097/ALN.0b013e3181fcf590
Liu
,
C.
,
Edwards
,
S.
,
Gong
,
Q.
,
Roberts
,
N.
, &
Blumhardt
,
L. D.
(
1999
).
Three dimensional MRI estimates of brain and spinal cord atrophy in multiple sclerosis
.
Journal of Neurology, Neurosurgery, and Psychiatry
,
66
(
3
),
323
330
.
Longstaff
,
A.
(
2000
).
Neuroscience.
New York
:
Springer
.
MacGregor
,
R. J.
(
1987
).
Neural and brain modeling.
San Diego, CA
:
Academic Press
.
MacGregor
,
R. J.
(
1993
).
Theoretical mechanics of biological neural networks.
Boston
:
Academic Press
.
Marc
,
R. E.
,
Jones
,
B. W.
,
Lauritzen
,
J. S.
,
Watt
,
C. B.
, &
Anderson
,
J. R.
(
2012
).
Building retinal connectomes
.
Current Opinions in Neurobiology
,
22
(
4
),
568
574
. doi:10.1016/j.conb.2012.03.011
Markram
,
H.
(
2006
).
The Blue Brain Project
.
Nature Reviews Neuroscience
,
7
(
2
),
153
160
. doi:10.1038/nrn1848
McLean
,
J. H.
,
Shipley
,
M. T.
,
Nickell
,
W. T.
,
Aston-Jones
,
G.
, &
Reyher
,
C. K.
(
1989
).
Chemoanatomical organization of the noradrenergic input from locus coeruleus to the olfactory bulb of the adult rat
.
Journal of Computational Neurology
,
285
(
3
),
339
349
. doi:10.1002/cne.902850305
Meinertzhagen
,
I. A.
(
2016
).
Connectome studies on drosophila: A short perspective on a tiny brain
.
Journal of Neurogenetics
,
30
(
2
),
62
68
. doi:10.3109/01677063.2016.1166224
Melzack
,
R.
, &
Wall
,
P. D.
(
1965
).
Pain mechanisms: A new theory
.
Science
,
150
(
3699
),
971
979
.
Modha
,
D. S.
, &
Singh
,
R.
(
2010
).
Network architecture of the long-distance pathways in the macaque brain
.
Proceedings of the National Academy of Sciences
,
107
(
30
),
13485
13490
. doi:10.1073/pnas.1008054107
Netter
,
F. H.
(
2014
).
Atlas of human anatomy.
Philadelphia, PA
:
Saunders/Elsevier
.
Noback
,
C. R.
(
2005
).
The human nervous system: Structure and function.
Totowa, NJ
:
Humana Press
.
Nógrádi
,
A.
, &
Vrbová
,
G.
(
2013
).
Anatomy and physiology of the spinal cord
, In
Madame Curie Bioscience Database [Internet].
Austin, TX
:
Landes Bioscience
.
Noori
,
H. R.
,
Schottler
,
J.
,
Ercsey-Ravasz
,
M.
,
Cosa-Linan
,
A.
,
Varga
,
M.
,
Toroczkai
,
Z.
, …
Spanagel
,
R.
(
2017
).
A multiscale cerebral neurochemical connectome of the rat brain
.
PLoS Biol.
,
15
(
7
),
e2002612
. doi:10.1371/journal.pbio.2002612
Park
,
H. G.
, &
Carmel
,
J. B.
(
2016
).
Selective manipulation of neural circuits
.
Neurotherapeutics
,
13
(
2
),
311
324
. doi:10.1007/s13311-016-0425-7
Piatkevich
,
K. D.
,
Jung
,
E. E.
,
Straub
,
C.
,
Linghu
,
C.
,
Park
,
D.
,
Suk
,
H. J.
, …
Boyden
,
E. S.
(
2018
).
A robotic multidimensional directed evolution approach applied to fluorescent voltage reporters
.
Nature Chemical Biology
,
14
(
4
),
352
360
. doi:10.1038/s41589-018-0004-9
Poldrack
,
R. A.
, &
Farah
,
M. J.
(
2015
).
Progress and challenges in probing the human brain
.
Nature
,
526
(
7573
),
371
379
. doi:10.1038/nature15692
Poldrack
,
R. A.
, &
Yarkoni
,
T.
(
2016
).
From brain maps to cognitive ontologies: Informatics and the search for mental structure
.
Annu. Rev. Psychol.
,
67
,
587
612
. doi:10.1146/annurev-psych-122414-033729
Ramirez-Jarquin
,
U. N.
, &
Tapia
,
R.
(
2018
).
Excitatory and inhibitory neuronal circuits in the spinal cord and their role in the control of motor neuron function and degeneration
.
ACS Chemical Neuroscience
. doi:10.1021/acschemneuro.7b00503
Ramirez-Mahaluf
,
J. P.
,
Roxin
,
A.
,
Mayberg
,
H. S.
, &
Compte
,
A.
(
2015
).
A computational model of major depression: The role of glutamate dysfunction on cingulo-frontal network dynamics
.
Cereb. Cortex
,
9
,
1
20
. doi:10.1093/cercor/bhv249
Remahl
,
S.
,
Cullheim
,
S.
, &
Ulfhake
,
B.
(
1985
).
Dimensions and branching patterns of triceps surae alpha-motor axons and their recurrent axon collaterals in the spinal cord during the postnatal development of the cat
.
Brain Res.
,
355
(
2
),
193
200
.
Roberts
,
A.
,
Conte
,
D.
,
Hull
,
M.
,
Merrison-Hort
,
R.
,
al Azad
,
A. K.
,
Buhl
,
E.
, …
Soffe
,
S. R.
(
2014
).
Can simple rules control development of a pioneer vertebrate neuronal network generating behavior?
J. Neurosci.
,
34
(
2
),
608
621
. doi:10.1523/JNEUROSCI.3248-13.2014
Ryan
,
K.
,
Lu
,
Z.
, &
Meinertzhagen
,
I. A.
(
2016
).
The CNS connectome of a tadpole larva of Ciona intestinalis (L.) highlights sidedness in the brain of a chordate sibling
.
Elife
,
5
. doi:10.7554/eLife.16962
Salkoff
,
D. B.
,
Zagha
,
E.
,
Yuzgec
,
O.
, &
McCormick
,
D. A.
(
2015
).
Synaptic mechanisms of tight spike synchrony at gamma frequency in cerebral cortex
.
J. Neurosci.
,
35
(
28
),
10236
10251
. doi:10.1523/JNEUROSCI.0828-15.2015
Sandberg
,
A.
, &
Bostrom
,
N.
(
2008
).
Whole brain emulation: A roadmap
. (
Technical Report 2008-3
).
Oxford
:
Oxford University
.
Sharma
,
A.
,
Wolf
,
D. H.
,
Ciric
,
R.
,
Kable
,
J. W.
,
Moore
,
T. M.
,
Vandekar
,
S. N.
, …
Satterthwarte
,
T. D.
(
2017
).
Common dimensional reward deficits across mood and psychotic disorders: A connectome-wide association study
.
American Journal of Psychiatry
,
174
(
7
),
657
666
. doi:10.1176/appi.ajp.2016.16070774
Shortland
,
P.
, &
Woolf
,
C. J.
(
1993
).
Morphology and somatotopy of the central arborizations of rapidly adapting glabrous skin afferents in the rat lumbar spinal cord
.
J. Comp. Neurol.
,
329
(
4
),
491
511
. doi:10.1002/cne.903290406
Shortland
,
P.
,
Woolf
,
C. J.
, &
Fitzgerald
,
M.
(
1989
).
Morphology and somatotopic organization of the central terminals of hindlimb hair follicle afferents in the rat lumbar spinal cord
.
J. Comp. Neurol.
,
289
(
3
),
416
433
. doi:10.1002/cne.902890307
Siekmeier
,
P. J.
, &
vanMaanen
,
D. P.
(
2013
).
Development of antipsychotic medications with novel mechanisms of action based on computational modeling of hippocampal neuropathology
.
PLoS ONE
,
8
(
3
),
e58607
. doi:10.1371/journal.pone.0058607
Silva
,
A. C.
(
2017
).
Anatomical and functional neuroimaging in awake, behaving marmosets
.
Developmental Neurobiology
,
77
(
3
),
373
389
. doi:10.1002/dneu.22456
Smith
,
K. M.
,
Boyle
,
K. A.
,
Mustapa
,
M.
,
Jobling
,
P.
,
Callister
,
R. J.
,
Hughes
,
D. I.
, …
Graham
,
B. A.
(
2016
).
Distinct forms of synaptic inhibition and neuromodulation regulate calretinin-positive neuron excitability in the spinal cord dorsal horn
.
Neuroscience
,
326
,
10
21
. doi:10.1016/j.neuroscience.2016.03.058
Stobb
,
M.
,
Peterson
,
J. M.
,
Mazzag
,
B.
, &
Gahtan
,
E.
(
2012
).
Graph theoretical model of a sensorimotor connectome in zebrafish
.
PLoS ONE
,
7
(
5
),
e37292
. doi:10.1371/journal.pone.0037292
Szigeti
,
B.
,
Gleeson
,
P.
,
Vella
,
M.
,
Khayrulin
,
S.
,
Palyanov
,
A.
,
Hokanson
,
J.
, …
Larson
,
S.
(
2014
).
Openworm: An open-science approach to modeling Caenorhabditis elegans
.
Frontiers in Computational Neuroscience
,
8
,
137
. doi:10.3389/fncom.2014.00137
Tang
,
Y.
,
Sun
,
W.
,
Toga
,
A. W.
,
Ringman
,
J. M.
, &
Shi
,
Y.
(
2017
).
A probabilistic atlas of human brainstem pathways based on connectome imaging data
.
Neuroimage
,
169
,
227
239
. doi:10.1016/j.neuroimage.2017.12.042
Tinaz
,
S.
,
Lauro
,
P. M.
,
Ghosh
,
P.
,
Lungu
,
C.
, &
Horovitz
,
S. G.
(
2017
).
Changes in functional organization and white matter integrity in the connectome in Parkinson's disease
.
Neuroimage: Clinical
,
13
,
395
404
. doi:10.1016/j.nicl.2016.12.019
Usoskin
,
D.
,
Furlan
,
A.
,
Islam
,
S.
,
Abdo
,
H.
,
Lonnerberg
,
P.
,
Lou
,
D.
, …
Ernfors
,
P.
(
2015
).
Unbiased classification of sensory neuron types by large-scale single-cell RNA sequencing
.
Nature Neuroscience
,
18
(
1
),
145
153
. doi:10.1038/nn.3881
van der Horn
,
H. J.
,
Kok
,
J. G.
,
de Koning
,
M. E.
,
Scheenen
,
M. E.
,
Leemans
,
A.
,
Spikman
,
J. M.
, …
van der Naalt
,
J.
(
2017
).
Altered wiring of the human structural connectome in adults with mild traumatic brain injury
.
Journal of Neurotrauma
,
34
(
5
),
1035
1044
. doi:10.1089/neu.2016.4659
Van Essen
,
D. C.
,
Smith
,
S. M.
,
Barch
,
D. M.
,
Behrens
,
T. E.
,
Yacoub
,
E.
,
Ugurbil
,
K.
, …
W. U-Minn HCP Consortium
. (
2013
).
The Wu-Minn Human Connectome Project: An overview
.
NeuroImage
,
80
,
62
79
. doi:10.1016/j.neuroimage.2013.05.041
Wang
,
Z.
,
Dai
,
Z.
,
Gong
,
G.
,
Zhou
,
C.
, &
He
,
Y.
(
2015
).
Understanding structural-functional relationships in the human brain: A large-scale network perspective
.
Neuroscientist
,
21
(
3
),
290
305
. doi:10.1177/1073858414537560
Wheeler-Kingshott
,
C. A.
,
Stroman
,
P. W.
,
Schwab
,
J. M.
,
Bacon
,
M.
,
Bosma
,
R.
,
Brooks
,
J.
, …
Tracey
,
I.
(
2014
).
The current state-of-the-art of spinal cord imaging: Applications
.
NeuroImage
,
84
,
1082
1093
. doi:10.1016/j.neuroimage.2013.07.014
Willis
,
W. D.
, &
Coggeshall
,
R. E.
(
2004a
).
Sensory mechanisms of the spinal cord: Ascending tracts and their descending control.
New York
:
Kluwer Academic/Plenum
.
Willis
,
W. D.
, &
Coggeshall
,
R. E.
(
2004b
).
Sensory mechanisms of the spinal cord: Primary afferent neurons and the spinal dorsal horn.
New York
:
Kluwer Academic/Plenum
.
Woolf
,
C. J.
(
1987
).
Central terminations of cutaneous mechanoreceptive afferents in the rat lumbar spinal cord
.
Journal of Comparative Neurology
,
261
(
1
),
105
119
. doi:10.1002/cne.902610109
Yap
,
C. B.
(
1967
).
Spinal segmental and long-loop reflexes on spinal motoneurone excitability in spasticity and rigidity
.
Brain
,
90
,
887
896
.
Zhang
,
T. C.
,
Janik
,
J. J.
,
Peters
,
R. V.
,
Chen
,
G.
,
Ji
,
R. R.
, &
Grill
,
W. M.
(
2015
).
Spinal sensory projection neuron responses to spinal cord stimulation are mediated by circuits beyond gate control
.
Journal of Neurophysiology
,
114
(
1
),
284
300
. doi:10.1152/jn.00147.2015