SIGNAL PROPAGATION VIA OPEN-LOOP INTRATHALAMIC ARCHITECTURES: A COMPUTATIONAL MODEL.

Propagation of signals across the cerebral cortex is a core component of many cognitive processes and is generally thought to be mediated by direct intracortical connectivity. The thalamus, by contrast, is considered to be devoid of internal connections and organized as a collection of parallel inputs to the cortex. Here, we provide evidence that "open-loop" intrathalamic pathways involving the thalamic reticular nucleus (TRN) can support propagation of oscillatory activity across the cortex. Recent studies support the existence of open-loop thalamo-reticulo-thalamic (TC-TRN-TC) synaptic motifs in addition to traditional closed-loop architectures. We hypothesized that open-loop structural modules, when connected in series, might underlie thalamic and, therefore cortical, signal propagation. Using a supercomputing platform to simulate thousands of permutations of a thalamocortical network based on physiological data collected in mice, rats, ferrets, and cats and in which select synapses were allowed to vary both by class and individually, we evaluated the relative capacities of closed- and open-loop TC-TRN-TC synaptic configurations to support both propagation and oscillation. We observed that 1) signal propagation was best supported in networks possessing strong open-loop TC-TRN-TC connectivity; 2) intrareticular synapses were neither primary substrates of propagation nor oscillation; and 3) heterogeneous synaptic networks supported more robust propagation of oscillation than their homogeneous counterparts. These findings suggest that open-loop, heterogeneous intrathalamic architectures might complement direct intracortical connectivity to facilitate cortical signal propagation.Significance Statement: Interactions between the dorsal thalamus and thalamic reticular nucleus (TRN) are speculated to contribute to phenomena such as arousal, attention, sleep, and seizures. Despite the importance of the TRN, the synaptic microarchitectures forming the basis for dorsal thalamus-TRN interactions are not fully understood. The computational neural model we present incorporates "open-loop" thalamo-reticular-thalamic (TC-TRN-TC) synaptic motifs, which have been experimentally observed. We elucidate how open-loop motifs possess the capacity to shape the propagative properties of signals intrinsic to the thalamus and evaluate the wave dynamics they support relative to closed-loop TC-TRN-TC pathways and intrareticular synaptic connections. Our model also generates predictions regarding how different spatial distributions of reticulothalamic and intrareticular synapses affect these signaling properties.


Introduction
Propagation of activity across the cerebral cortex is thought to underlie multiple cognitive processes, as well as pathologic processes such as epilepsy and migraine (Leao, 1944;Muller et al., 2014Muller et al., , 2016Kokkinos et al., 2017). Cortical regions are highly interconnected via direct axonal projections as well as via polysynaptic pathways involving the basal ganglia and thalamus (Parent and Hazrati, 1995;Theyel et al., 2010). Cortical signal propagation is generally thought to be mediated via direct cortical connections (Felleman and Van Essen, 1991;Kötter and Sommer, 2000), but recent evidence suggests that the thalamus serves as a control point to modify cortical activity during cognitive processes such as attentional shifting (Wimmer et al., 2015). An advantage of a thalamic mode of signal propagation is the efficiency by which modulatory influences may control thalamic, and therefore cortical, propagation. The thalamus, however, is generally thought to have limited internal connectivity and therefore limited capacity to serve as a substrate for signal propagation.
Based on previous studies of open-loop TC-TRN-TC synaptic organization, we hypothesized that open-loop synaptic modules might underlie intrathalamic and therefore intracortical signal propagation. Accordingly, we systematically examined thousands of permutations of a novel network model comprising thalamic, reticular, and thalamorecipient, layer-4 cortical (Co) neurons to evaluate the efficacy of known thalamic synaptic motifs (open-loop pathways, closed-loop pathways, and chemical and electrical intrareticular synapses), in isolation and in combination, in mediating signal transmission across the thalamus and cortex.

Network architecture and simulations
We constructed a baseline model network based on Willis et al. (2015) by connecting in series three thalamocortical pathways, each consisting of a TC, TRN, and layer-4 Co neuron (for a three-by-three neuron network); the physiological data used in our model were collected in mice, rats, ferrets, and cats of both sexes. Permutations of the baseline network potentially featured both closed-loop and/or open-loop TC-TRN-TC motifs, with the latter constituting one mode of connectivity between parallel thalamocortical pathways. Intrareticular synapses represented the other connections between pathways, based on the identification of both GABAergic (Ahlsén and Lindström, 1982;Steriade et al., 1990;Cox et al., 1996;Sanchez-Vives et al., 1997;Shu and McCormick, 2002;Deleuze and Huguenard, 2006;Lam et al., 2006) and electrical synapses (Landisman et al., 2002;Fuentealba et al., 2004;Long et al., 2004;Deleuze and Huguenard, 2006;Lam et al., 2006) between TRN neurons. Thus, we included three different polysynaptic configurations between thalamocortical pathways in our network (Fig. 1B, from left to right): (1) those with a GABAergic intrareticular synapse (TRN-TRN GABA ), (2) those with an electrical intrareticular synapse (TRN-TRN Elec ), and (3) open-loop TC-TRN-TC pathways. Thalamic, reticular, and cortical cell layers were aligned topographically, such that TC 1 projected to both TRN 1 and Co 1 (Jones, 1975;Steriade et al., 1993;Destexhe et al., 1998;Sohal et al., 2000;Sherman and Guillery, 2001). The divergence of thalamic and reticular synapses in the model was constrained to accommodate open-loop TC-TRN-TC architectures, which depend on a lack of recurrent feedback to the downstream TC neuron, in a subset of the simulated network variants: to this end, every TC neuron projected to exactly one TRN neuron, while single TRN neurons could project to either one TC neuron (whether recurrently or laterally), as in the case of entirely closed-loop or open-looped TC-TRN-TC motifs, or two TC neurons (one recurrently and one laterally), if participating in a pathway expressing some intermediate degree of openness (Fig. 1C).
To analyze how each variety of interpathway connection contributed to network dynamics, permutations of the baseline network were generated by varying three synaptic properties associated with each of the interpathway synaptic motifs; moreover, these parameters were either varied in a homogeneous or heterogeneous manner. In the case of homogeneously varied synaptic network permutations, the synaptic parameters associated with three interpathway motifs varied uniformly as a class, with all external, TC-TRN, and thalamocortical (TC-Co) synaptic conductances held constant: (1) TRN-TRN GABA synapses ranged in conductance between 0 and 450 nS, (2) TRN-TRN Elec synapses ranged in coupling coefficient between 0 and 0.36, and (3) a TC-TRN-TC "openness" Three possible pathways through which a signal might propagate from one TC neuron to another via the TRN. C, Baseline thalamocortical model network. Broken-line synapses were allowed to vary either as a class (homogeneously) or independently of one another (heterogeneously). The black arrow corresponds to the fixed, external stimulus applied to TC 1 . D, Sample cortical spike histograms (detrended) in a network permutation responding to a fixed, sustained stimulus delivered to TC 1 (black bar beneath the lowest trace). The propagation score assigned to any network permutation was quantified as the amplitude of the initial stimulus-evoked response in the detrended Co 3 histogram; response propagation across the cortical subnetwork (orange arrow) was consistently linear, and thus the initial response in Co 3 was observed at a fixed interval relative to the onset of stimulation. Oscillation intrinsic to any network variant was quantified as the amplitude of the first off-center peak in the normalized autocorrelogram (right) of poststimulation activity in the detrended Co 3 histogram (within the broken black box). The initial 400 ms of activity preceding the fixed stimulus (in gray) is shown here for each histogram but was not included in the calculations of either propagation or oscillation. Note that the bin heights in the Co 1 histogram shown here were truncated in order to maintain identical vertical scaling across all three cortical histograms. coefficient, defined as the weight distribution of lateral versus recurrent reticulothalamic (TRN-TC) connectivity, varied between 0 (completely closed-loop) and 1.0 (completely open-loop) and with a baseline TRN-TC conductance of 80 nS. Thus, for a network variant possessing an openness coefficient of 0.4 (i.e., exhibiting slightly more closed-loop than open-loop TC-TRN-TC connectivity), the laterally inhibitory TRN-TC synapses in the network, TRN 1 !TC 2 and TRN 2 !TC 3 , would carry a conductance of 32 nS (0.4 Â 80 nS), while the recurrent TRN-TC synapses, TRN 1 !TC 1 , TRN 2 !TC 2 , and TRN 3 !TC 3 , would exhibit a conductance of 48 nS (0.6 Â 80 nS). For the heterogeneously varied synaptic network variants, all TRN-TRN and TRN-TC synapses were allowed to vary independently. Domains for each of the synaptic variables were selected to include the range of conductance or coupling strengths reported in physiological measurements and/or used in similar neural models (Destexhe et al., 1996a(Destexhe et al., , 1998Sohal and Huguenard, 1998;Sohal et al., 2000;Landisman et al., 2002;Long et al., 2004;Traub et al., 2005).
Ongoing afferent synaptic input was delivered to every TC neuron in the model as Poisson-modulated spike trains centered at 40 Hz. An additional 200-Hz pulse train was applied to neuron TC 1 between t = 0.400 s and t = 1.500 s during every network simulation run. This high-frequency stimulus was modeled on those used to elicit spindle-like waves in a ferret thalamoreticular slice preparation Kim et al., 1995). A given network's output was compiled by assembling spike histograms (10-ms bins) averaging 1000 simulations for every cortical neuron (Fig.  1D). We quantified network dynamics as a function of variable TC-TRN-TC and intrareticular synaptic architectures by defining and measuring two properties inherent to stimulus-evoked responses in each network variant: propagation and oscillation, with the latter included in light of the fact that many characterized thalamic waveforms both oscillate and propagate through the thalamus and cortex (Sherman and Guillery, 2001). Network properties were quantified in the most downstream element of the cortical output layer, Co 3 . Propagation across a network was quantified as the amplitude of the initial stimulus-evoked response in the detrended Co 3 histogram. The degree of oscillation supported by each network permutation was defined as the amplitude of the first off-center peak in the normalized autocorrelogram of poststimulation activity (Fig. 1D). Both propagation and oscillation scores are reported as normalized to the maximum scores tabulated for each property. Given the high prevalence of propagating oscillatory waves in the cerebral cortex (for review, see Muller et al., 2018), we furthermore defined a composite "optimization" (Op) metric to measure the capacity of networks to simultaneously support and balance between propagation (Pr) and oscillation (Os): Op ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Pr 2 1Os 2 p À jPr À Osj: Intrinsic neuronal models Single-compartment TC, TRN, and cortical model neurons obeyed Hodgkin-Huxley kinetics, with membrane potentials V varying according to the first-order differential equation: where C is the membrane capacitance, g L and E L are the leakage conductance and reversal potential, respectively, and g i (V) and E i are the dynamic conductance and reversal potential, respectively, of the ith voltage-gated, ligandgated (chemical synaptic), or electrical synaptic conductance (for electrical synaptic conductances, the effective reversal potential is equal to the presynaptic membrane potential; see Eq. 3a). All three varieties of model neurons expressed both the standard transient sodium (I Na ) and delayed-rectifier potassium (I K ) currents. TC and TRN neurons additionally included a T-type calcium conductance (T-current; I T ) and hyperpolarization-activated cation current (H-current; I H ), following the TC model of Deleuze et al. (2012). Both TRN and layer-4 Co cells expressed a slow, non-inactivating potassium conductance (I M ), following the modeling of Pospischil et al. (2008), which accounts for the spike-frequency adaptation previously reported in physiological recordings from these neurons (Yamada et al., 1989;Willis et al., 2015). A list of intrinsic model cell parameters, including current conductances, reversal potentials, selected gating kinetics, and membrane capacitance, can be found in Table 1.

Synaptic models
The kinetics of chemical synapses in our model network conformed to the synaptic depression model of Tsodyks and Markram (1997). This model presupposes a finite quantity of "resources," akin to synaptic vesicles, capable of being released by the presynaptic neuron; these resources can exist in an active, inactive, or recovered state. A parameter U SE characterizes the fraction of recovered resources that can be converted to an active state (i.e., for release by the presynaptic neuron) following action potential induction in the presynaptic axon terminal (s). Following resource activation, synapses inactivate according to the time constant t inact ; resources become available again for activation after a recovery period described by the time constant t recov . These parameters, along with the neurotransmitters, postsynaptic conductances, and reversal potentials characterizing all of the chemical synapses in our model, are given in Table 2.
Glutamatergic TC-TRN and TC-Co and baseline GABAergic TRN-TC synaptic parameters matched those of Willis et al. (2015), with the latter synapses allowed to vary in conductance as described above. TRN-TC signaling was mediated exclusively through GABA A receptors, mirroring other thalamic and thalamocortical models in which the slower TRN-TC GABA B conductance was omitted (Traub et al., 2005;Rogala et al., 2013;Pham and Haas, 2018). Although evidence has been presented challenging the existence of GABAergic intrareticular synapses in certain mammalian species and age groups (Pinault et al., 1997;Landisman et al., 2002;Pinault, 2004;Cruikshank et al., 2010;Hou et al., 2016), our model avoided making assumptions regarding their presence, strength, or spatial distribution by allowing the associated synaptic conductances to vary over a range of physiological values, including zero, and in distribution. The reversal potential, conductance, and kinetics of the external synapses projecting to the TC neurons were directly based on retinogeniculate synapses (Chen and Regehr, 2003), although the generic nature of the external inputs in our model allows them to represent not only immediately upstream sensory input but also brainstem modulation (e.g., serotonergic, adrenergic) known to act on thalamic nuclei (Siegel and Sapru, 2015).
Electrical synapses between TRN neurons were based on the Cx36-dependent intrareticular gap junctions first identified by Landisman et al. (2002). For TRN neurons, the sum of electrical synaptic currents (I Elec ) entering any postsynaptic neuron j from presynaptic neuron(s) i was included in the rightmost term from Equation 2 and calculated as where g ij , the gap junction conductance, was itself calculated as where CC was the electrical coupling coefficient between TRN neurons i and j, g j was the membrane conductance of the postsynaptic neuron, and D(x) was a scaling factor that depended on the physical distance between the coupled TRN neurons (Dayan and Abbott, 2005;Shimizu and Stopfer, 2013). Although the individual gap junctions comprising TRN-TRN Elec synapses used in the model were not explicitly coded for, differences in TRN-TRN coupling between different electrical synapses and network permutations implicitly reflected differing gap junction densities, following Traub et al. (2005). TRN-TRN Elec synapses were symmetrical (non-rectifying), such that g ij =g ji . We extrapolated the attenuation of intrareticular synaptic strength as a function of intracellular distance based on mappings of intrinsic connections within the TRN along a horizontal (anteroposterior) plane assembled by Deleuze and Huguenard (2006). Assuming (1) an intracellular distance of 50 mm between adjacent TRN neurons, (2) a distance x (in multiples of 50 mm) between non-adjacent neurons, and (3) a Gaussian falloff in synaptic strength (Sohal et al., 2000), the baseline (adjacent-neuron) conductances of TRN-TRN GABA and TRN-TRN Elec  synapses were scaled for non-adjacent synapses using the function where l GABA = 531 mm and l Elec = 130 mm. Given the small spatial scale of our model, synaptic delays associated with finite axonal conductance times within the TRN and between the TRN and dorsal thalamus were disregarded, mirroring the simplification incorporated into previous thalamic and thalamocortical models simulating synaptic interactions on the order of 100 mm (Golomb et al., 1996;Traub et al., 2005). Although small (;1 ms) thalamocortical delays were inserted into the network model of Traub et al. (2005), these were likewise omitted on the basis of the cortex functioning solely as an output layer in our model.

Quantification and statistical analysis
Our model was coded, simulated, and analyzed in MATLAB R2018b (MathWorks) using both a Dell Inspiron 3847 and Hewlett-Packard Z840 running Windows 10 and nodes on the Illinois Campus Cluster (National Center for Supercomputing Applications, University of Illinois at Urbana-Champaign). Simulations, of which there were 1000 for every network permutation, employed 0.1-ms time steps, with temporal integration based on the hybrid analytic-numeral integration method of Moore and Ramon (1974), which optimizes between accurate solutions to Hodgkin-Huxley and synaptic models and computational efficiency. All simulations commenced with a 200-ms equilibration period, during which no external stimulation was delivered to TC neurons; this allowed all network elements to attain steady-state conditions. The number of homogeneously and heterogeneously varied synaptic network variants generated were 770 and 12,681, respectively.
Statistical analysis was performed in both MATLAB and R (R Core Team, 2013), with the glmnet package (Friedman et al., 2010) used within the latter platform to perform regression analyses. Multiple linear regression was employed to establish rudimentary relationships between synaptic classes (homogeneously synaptic networks) or individual synapses (heterogeneously synaptic networks) and each of the two studied network properties, even in instances where these relationships deviated from linearity. Second-order polynomial (2°) regression models with interaction terms elucidated how synaptic interactions and nonlinearities affected these network properties. Regressions were optimized using elastic net regularization, with the specific regularization hyperparameter a selected to minimize each regression model's root-mean-square error (RMSE). To convey the relative influence of different synaptic classes or individual synapses on dynamic network properties, all regression coefficients are reported here as normalized to the coefficient with the largest absolute value; the effects corresponding to normalized regression coefficients (NRCs) with absolute values of ,0.05 were disregarded as negligibly influential on network dynamics. Both unpaired Student's t tests and one-way analysis of variance models were used to compare the mean property scores between different sets of networks, with Tukey's honestly significant difference tests used to ascertain pairwise difference between groups in the latter; standard errors of the mean (SEMs) were used as a measure of variance, and null hypotheses were rejected at probability values (p values) below 0.05. Kolmogorov-Smirnov and Levene's tests were employed to confirm normality and homogeneity of variance, respectively, when using parametric meancomparison tests; data were transformed as needed to conform to these prerequisites.

Code accessibility
The code/software described in the paper is freely available online at https://github.com/JeffreyWBrown/ Open-loop-TC-TRN-TC. The code is available as Extended Data 1.

Propagation and oscillation in homogeneously varied synaptic models
Stimulus-evoked responses propagated linearly across the length of homogeneous synaptic networks, occurring at average fixed intervals of 93.31 6 0.35 ms (mean 6 SEM; range, 60-110 ms) between adjacent thalamocortical pathways, across all model permutations and with a mean velocity of 0.54 mm/s, assuming a 50-mm separation between adjacent neurons in each network layer. All 770 homogeneous network variants were ranked according to their cortical propagation scores ( Fig. 2A, top). TC neurons exhibited both tonic firing and bursting activity, with the former mode more frequently observed (Fig. 2B).
Multiple linear regression analysis (R 2 = 0.793, RMSE = 0.047, p , 0.0001) demonstrated a strong positive correlation between the TC-TRN-TC openness coefficient and propagation score (NRC = 1.000). By contrast, chemical and electrical TRN-TRN synaptic connectivity tended to modestly diminish propagation (NRC = À0.173 and NRC = À0.136, respectively; Table 3). Further, other excitatory connectivity, such as corticocortical or corticothalamic connectivity, often postulated as being important for cortical signal propagation (Felleman and Van Essen, 1991;Kötter and Sommer, 2000;Theyel et al., 2010), was not necessary. Thus, the homogeneously varied synaptic network permutations that best accommodated signal propagation were generally ones with weak or absent synapses between TRN neurons and strong open-loop TC-TRN-TC connections. For example, Network a, which epitomizes this architecture, exhibited robust signal propagation in response to a fixed stimulus delivered to TC 1 ; a representative simulation of this network is shown in Figure 2B, left, and its position within the parameter space depicted in the left-sided heat map of Figure 2C is labeled. Stimulus-evoked activity in this network tended to propagate efficiently from Co 1 to Co 3 : near-synchronous propagation cascades were elicited in both the TRN and cortical layers of the model, having been stimulated by propagating activity in upstream TC neurons. Smooth,   linear propagation of action potentials across the network depended on the synchronous induction of IPSPs and the ensuing postinhibitory rebound spikes in TC neurons, as mediated by T-type Ca 21 channels and driven by inhibition from the TRN, which occurred reliably and at fixed intervals in Network a. Relative to Network a, other network permutations exhibiting stronger intrareticular synapses did not support propagation as efficiently. We surmise that TRN-TRN GABA synaptic connections reduced the incidence of IPSPs in TC neurons required for signal propagation across the network, while electrical coupling between TRN neurons destructively shunted a propagating signal away from the thalamoreticular lattice through which it predominantly traversed the network. A 2°regression model of propagation as a function of all three synaptic class variables (R 2 = 0.842, RMSE = 0.041, p , 0.0001; Table 3) revealed a modestly negative interaction term between TRN-TRN Elec synapses and TC-TRN-TC openness (NRC = À0.365), indicating that in networks where both electrical synapses were strong and TC-TRN-TC openness high, the extent of supported propagation diminished nonlinearly; a smaller negative interaction between TRN-TRN GABA synapses and TC-TRN-TC openness was also observed (NRC = À0.152). Together, these terms suggested that propagation was more significantly affected by connections in the TRN layer within the more open-loop networks. This relationship was evident, for example, in the right-sided heat map of Figure 2C, in which propagation scores more markedly decreased with increasing TRN-TRN electrical coupling as TC-TRN-TC openness itself increased.
Oscillatory responses recurred in Co 3 neurons at a mean frequency of 9.07 6 0.2 Hz (range, 7.14-12.50 Hz) across all homogeneous model permutations. Propagation and oscillation scores across all 770 homogeneous networks were strongly anticorrelated (Pearson's r = À0.671, p , 0.0001). Accordingly, oscillation was best accommodated in network permutations exhibiting strongly closed-loop connectivity ( Fig. 2A, bottom); however, the capacity to support oscillation was neither markedly linear nor monotonically decreasing as a function of increasing openness coefficient (Fig. 2D). Rather, a one-way analysis of variance with Tukey's tests revealed that, on average, oscillation scores peaked and remained statistically indistinguishable from one another across the subset of network permutations with openness coefficients between 0 and 0.4, with scores then decreasing in a roughly linear fashion with increasing TC-TRN-TC openness (F (10,759) = 137.8, p , 0.0001); this result was consistent with a 2°regression model of oscillation (R 2 = 0.630, RMSE = 0.128, p , 0.0001; Table 3), in which the linear and quadratic terms in TC-TRN-TC openness were associated with NRCs of À1.000 and À0.052, respectively, and the effects of TRN-TRN GABA (NRC = 0.060) and TRN-TRN Elec (NRC regularized to 0) on oscillation were weakly positive and negligible, respectively. Taken with the analysis of propagation, these data suggest that networks with mixed open-and closed-loop connectivity (which is likely close to physiological reality) can support the coexistence of oscillation and propagation (see below, Propagation and oscillation in heterogeneously varied synaptic models).
The predominant mechanism by which oscillation arose in Co 3 was through postinhibitory rebound in TC 3 , as engendered by the strong recurrent inhibition found in network permutations exhibiting primarily closed-loop TC-TRN-TC connectivity. This mode of oscillation was exemplified by Network b , a strongly closed-loop network variant. In the simulation shown of this network (Fig.  2B, right), oscillatory activity was enabled by a single epoch of signal propagation. Notably, neither the presence of strong GABAergic nor electrical intrareticular synapses in Network b exerted much effect on its ability to support oscillation, as predicted by the 2°regression model.

Propagation and oscillation in heterogeneously varied synaptic models
Recent studies have highlighted heterogeneity in TRN neuronal connectivity, synaptic physiology and chemical identities (Lee et al., 2007;Halassa et al., 2014;Clemente-Perez et al., 2017). We therefore examined the impact of allowing all synaptic connections involving the TRN to be independently varied. We constructed circuit-level schematics of linear regression models for propagation (Fig.  3A, top) and oscillation (Fig. 3A, bottom) as functions of the 14 synaptic variables in heterogeneous networks.
Propagation in heterogeneously varied synaptic networks increased chiefly as a function of increasing the strength of the more downstream of the two laterally inhibitory TRN-TC synapses, TRN 2 !TC 3 : the corresponding term in a linear regression model of propagation (R 2 = 0.742, RMSE = 0.069, p , 0.0001; Table 4) possessed an NRC of 1.000 (Fig. 3A, top). Propagation scores also scaled to a lesser extent with the more upstream laterally inhibitory reticulothalamic synapse, TRN 1 !TC 2 (NRC = 0.608). The two inhibitory intrareticular synapses originating at the rightmost end of the model network, TRN 3 !TRN 1 and TRN 3 !TRN 2 , both exerted a small negative effect on propagation (NRC = À0.087 and NRC = À0.084, respectively). Additionally, two TRN-TRN Elec synapses, TRN 1 =TRN 2 and TRN 1 =TRN 3 (where the "=" denotes an electrical synapses), marginally decremented propagation in heterogeneous networks, with NRCs of À0.051 and À0.072, respectively. These findings at an individual synaptic level comported with the observation that strong TRN-TRN interactions, whether chemical or electrical, tended to impede signal propagation in homogeneous network variants.
A 2°regression model (R 2 = 0.857, RMSE = 0.051, p , 0.0001; Table 4) disclosed a large, propagation-enhancing interaction between the two laterally inhibitory synapses (NRC = 0.753), underscoring the same dependence of propagation on strong open-loop TC-TRN-TC connectivity as seen in homogeneously synaptic networks, but additionally demonstrating that propagation scores increased nonlinearly as a function of simultaneously increasing the weights of TRN 1 !TC 2 and TRN 2 !TC 3 . Interactions between TRN-TRN synapses of either variety and TRN-TC synapses tended diminish propagation, as did those between recurrent and lateral inhibitory TRN-TC synapses. Taken together, the linear and 2°regression models indicated that heterogeneous network permutations with strong laterally inhibitory TRN-TC synapses tended to best support propagation. Consistent response propagation across the length of the network was epitomized by Network a', in which TRN 1 !TC 2 and TRN 2 !TC 3 were both relatively strong and those synapses impeding propagation comparatively weak (Fig. 3B, left).
We analyzed the relative capacities of homogeneously and heterogeneously varied synaptic networks to support propagation, oscillation, and optimization by comparing the 20 highest scores achieved by homogeneous and heterogeneous network permutations with respect to each performance metric. No significant differences in mean propagation scores between top-performing homogeneous and heterogeneous networks were disclosed (unpaired t test, t (38) = 0.46, p = 0.647; Fig. 4). We attributed this lack of differences to the fact that network permutations in which the synapses TRN 1 !TC 2 and TRN 2 !TC 3 were both maximally weighted would be equally capable of supporting robust signal propagation, regardless of whether these synapses were varied homogeneously or heterogeneously. By contrast, top-scoring heterogeneous network variants better supported both oscillation (t (38) = 13.88, p , 0.0001) and optimization (t (38) = 18.04, p , 0.0001) than their homogeneous counterparts. Because networks supporting the propagation of oscillatory activity would, by definition, score high with respect to optimization, these results not only confirmed that heterogeneous networks were more likely than homogeneous  Table 4. B, Representative simulations for two selected heterogeneous networks, whose normalized synaptic weights are depicted in the circuit diagrams. Networks a' and b ' respectively illustrate propagation and propagation of oscillation across the network. networks to accommodate this oscillatory mechanism, but furthermore disclosed that propagation of oscillation across the thalamocortical network was associated with higher oscillation scores than postinhibitory-driven oscillation in TC 3 , the predominant form of oscillation observed in homogeneous networks.

Discussion
The presented analysis suggests that open-loop TC-TRN-TC synaptic motifs (Fig. 1B, right) can function as a substrate for signal propagation across the thalamus, and by extension, cortical networks without the need for direct corticocortical, intrareticular or corticothalamic connectivity. Postinhibitory rebound mediated by T-type Ca 21 channels served as a mediator of both propagation and oscillation in the simulated networks. TRN-TRN connections, either chemical or electrical (Fig. 1B, left and middle), diminished horizontal propagation by disrupting the precise timing relationships required to propagate a signal across the network. Models with heterogeneously varied synapses outperformed those whose synapses varied as a class with respect to the propagation of oscillatory activity, consistent with the emerging literature documenting cellular and synaptic heterogeneity in the TRN (Lee et al., 2007;Halassa et al., 2014;Clemente-Perez et al., 2017). These data suggest that widespread propagating thalamic or thalamocortical activity, under both pathologic and physiological conditions, may be mediated, at least in part, by TC-TRN-TC connections. The model makes strong predictions that can be tested physiologically.
Additionally, the present model omitted explicit corticothalamic and corticoreticular synapses, both of which have been identified and physiologically characterized to varying degrees (Steriade et al., 1972;White and Hersch, 1982;De Curtis et al., 1989;Blumenfeld and McCormick, 2000;Zhang and Jones, 2004;Crandall et al., 2015), although the former were effectively amalgamated with both feedforward sensory and modulatory projections to the thalamus in the form of the Poisson-modulated external input we delivered to individual TC neurons. Both forms of feedback have been implicated in the spread of spindle waves and in the maintenance of their synchronization over large distance scales (on the order of the length of the mammalian forebrain) and are furthermore known to drive spindle wave formation and propagation in vivo by polysynaptically recruiting TC neurons via TRN-mediated postinhibitory rebound (Steriade et al., 1972;Roy et al., 1984;Suga and Ma, 2003;Sillito et al., 2006;Crandall et al., 2015;Sorokin et al., 2017). It should be noted, however, that short-range propagation of spindle waves, which can be elicited in isolated thalamic slice preparations Kim et al., 1995), is preserved following decortication, both in vivo and in silico Destexhe et al., 1998).
The signaling dynamics observed in our small-scale, broadly feedforward model would undoubtedly be altered by introducing corticothalamic and corticoreticular feedback, as well as corticocortical synapses: in particular, we can predict, based on other modeling studies that have systematically explored the contributions of such connections within fundamentally closed-loop thalamoreticular frameworks (Destexhe et al., 1998;Rogala et al., 2013), that cortical feedback to the dorsal thalamus and TRN would increase the frequency of postinhibitory rebound in relay cells of the former, while also increasing oscillatory activity through the introduction of new recurrent pathways between the thalamus and cortex. While the enhancement of postinhibitory rebound in the thalamus would a priori suggest an enhancement in the extent of signal propagation mediated through open TC-TRN-TC loops, our results demonstrate that the efficiency of such propagation can be attenuated by decohering activity introduced via intrareticular synapses: as such, whether corticothalamic, corticoreticular, and corticocortical synapses dynamics would ultimately enhance open-loopmediated propagation or decrement it by interfering with the temporal dynamics underlying this form of intrathalamic signaling would heavily depend on, among other things, the degree to which these additional synapses spatially diverged within the dorsal thalamus and TRN. Notwithstanding elaborations of cortical projections within our network model, our present results suggest that open-loop TC-TRN-TC architectures may mediate a novel form of intrinsic thalamic, and by extension, cortical signal propagation that exists independently of topdown modulation (for example, in localized regions of the thalamus where cortical innervation is sparse) and potentially in parallel to the modes of thalamocortical propagation in which both corticofugal pathways and corticocortical synapses are known to participate. Future efforts within our modeling paradigm stand to incorporate reciprocal cortical projections involving one or multiple cortical layers (Destexhe et al., 1998;Bazhenov et al., 2002;Traub et al., 2005;Izhikevich and Edelman, 2008;Rogala et al., 2013).

Comparison to related computational models and physiological data
Although the production of spindle waves was not an explicit objective of our study, some of the wave dynamics arising in our networks were nevertheless consistent with those inherent to spindle or spindle-like waves. Despite possessing higher degrees of TC!TRN and TRN!TC synaptic divergence and lacking the exclusively open-loop TC-TRN-TC architecture characterizing a subset of our network variants, other isolated thalamic models allowing for longitudinal wave propagation similarly accommodated this propagation along the lattice of interconnected TC and TRN neurons by way of laterally inhibitory TRN-TC synapses (Destexhe et al., 1996a;Golomb et al., 1996;Bazhenov et al., 1998); at short ranges, this mechanism of signal propagation also prevailed in larger-scale thalamocortical models, while corticothalamic projections acted to propagate activity to more distal sites (Destexhe et al., 1998;Destexhe and Sejnowski, 2003). Comparably, recurrently inhibitory TRN-TC synapses have been documented to play a vital role in the generation of oscillatory behavior in the thalamus (von Krosigk et al., 1993;Destexhe and Sejnowski, 2003). The temporal parameters of propagating and oscillation signals in our model also matched some of those previously reported: the mean signal propagation velocity and oscillation frequency measured across homogeneous networks fell within the ranges of spindle wave propagation velocities and intraspindle spike frequencies reported in both physiological and computational spindle wave studies (Andersen and Andersson, 1968;Steriade and Deschenes, 1984;Kim et al., 1995;Destexhe et al., 1996a;Golomb et al., 1996).
The TC neurons in our model network exhibited both tonic and bursting modes of firing, consistent with extensive physiological characterization (Sherman, 2001). The form of signal propagation generated in our networks via open-loop TC-TRN-TC synaptic pathways, which necessarily depended on TRN-driven postinhibitory rebound in the downstream TC neuron, could, in practice, be elicited regardless of whether the upstream TC neuron fired tonically or in bursts. However, bursting in thalamic neurons receiving open-loop inhibition from the TRN is associated with a greater fidelity of transmission to the cortex relative to tonic firing, as was systematically demonstrated in the computational model that directly inspired our present study (Willis et al., 2015); this finding furthermore holds relative to thalamocortical signaling efficiency more generally Reinagel et al., 1999, Swadlow andGusev, 2001;Krahe and Gabbiani, 2004).
One particularly notable point of departure relative to similar network models was the extent to which thalamoreticular, reticulothalamic, and thalamocortical synapses diverged. Although all three classes of synapses are known to diverge significantly and have been observed to target neuronal somata hundreds of microns from their origins (Jones, 1985;Cox et al., 1996Cox et al., , 1997Crabtree, 1996;Pinault and Deschênes, 1998;Alonso et al., 2001;Miller et al., 2001;Sherman and Guillery, 2001), the TC-TRN, TRN-TC, and TC-Co synapses in our model were constrained to remain strictly local and minimally divergent (or non-divergent, in the case of TC-TRN and TC-Co synapses). With respect to the first two classes of synapses, this constraint was imposed to probe the impact the disynaptic TC-TRN-TC open-loop motifs characterizing a subset of network permutations, which constituted one of the foci of our study, and analyze the signal propagation they may support. This neuroanatomical scheme contrasted with previous computational models featuring parallel, interconnected thalamoreticular pathways, in which both TC and TRN synapsed bidirectionally with several neighboring TRN and TC cells, respectively, within a radius of several hundred microns (Destexhe et al., 1996a(Destexhe et al., , 1998Golomb et al., 1996;Bazhenov et al., 1998;Sohal and Huguenard, 1998;Sohal et al., 2000;Traub et al., 2005;Izhikevich and Edelman, 2008). Furthermore, the limited synaptic divergence constrained the spatial and temporal scales over which propagating and oscillating signals persisted in our model. This does not necessarily imply that thalamic and/or cortical signal propagation mediated through open-loop TC-TRN-TC architectures would be inherently limited in either distance or duration, particularly when accounting for the comparative diversity and complexity in the spatial and temporal profiles of real sensory information integrated by the thalamus relative to the highly focal, time-fixed external stimulus approximations we employed to initialize responses reliably across simulations in our model networks, consistent with similarly simplified stimulus representations used in other thalamic or thalamocortical models (Destexhe et al., 1996a;Golomb et al., 1996;Bazhenov et al., 1998;Sohal and Huguenard, 1998;Traub et al., 2005). In light of the limits on the spatiotemporal coherence of signals intrinsic to our present model, however, we would not predict any qualitative changes in propagative or oscillatory dynamics were we to increase the length of our baseline network as presently constituted by adding in parallel additional TC, TRN, and Co neurons.

The functional implications of open-loop TC-TRN-TC synaptic motifs
The spread of activity from one cortical region to another is a foundational concept at the core of our understanding of sensory processing, higher-order cognitive functions such as attention and language, sleep-related oscillatory phenomena, and pathologic findings such as propagation of ictal discharges and migraine. It has long been speculated that the TRN could serve as a control point for large-scale cortical signal processing given its central location, the high degree of convergence of projections involved in attention, arousal, and emotion onto the TRN, and the TRN's particularly strong control over TC firing properties (Yingling and Skinner, 1976;Crick, 1984;Guillery et al., 1998;Brunia and Van Boxtel, 2001;McAlonan et al., 2006;John et al., 2013). Here, we showed that open-loop TC-TRN-TC architectures can support at least short-range thalamocortical signal propagation. Within the thalamus, these configurations have thus far been observed both within and across individual thalamic nuclei and are thought to serve as pathways for intra-and cross-modal modulation, respectively (Crabtree et al., 1998;Pinault and Deschênes, 1998;Crabtree and Isaac, 2002;Sherman, 2005, 2015;Kimura et al., 2007;Lee et al., 2010;Kimura, 2014); as has been previously surmised, these synaptic pathways could also plausibly lend themselves to sensory enhancement, multisensory integration, and attentional mechanisms (Crabtree and Isaac, 2002;Pinault, 2004;Willis et al., 2015;Crabtree, 2018). At a minimum, and as inferred from physiological studies, open-loop pathways should be fully capable of supporting signaling propagation from one thalamic relay neuron to another through a limited number of intervening synapses (with a disynaptic pathway serving as the shortest such configuration).
It should be emphasized that the specific functional roles of open-loop TC-TRN-TC pathways are likely to depend on their densities and distributions within the thalamus (Halassa and Acsády, 2016). If the morphologic, intrinsic, and synaptic heterogeneity of TRN neurons are any indication (Scheibel and Scheibel, 1966;Jones, 1975;Spreafico et al., 1991;Cox et al., 1996;Lee et al., 2007;Halassa et al., 2014;Clemente-Perez et al., 2017), it is reasonable to assume that both TC-TRN and TRN-TC synapses are distributed in a broadly heterogeneous manner across the thalamus. As underscored by our analysis, such synaptic heterogeneity is seemingly a prerequisite for the propagation of oscillatory signals, which, in the case of spindling, can occur in the thalamus independently of cortical involvement (von Krosigk et al., 1993;Bal et al., 1995;Kim et al., 1995) and necessarily involves both recurrently and laterally projecting TRN-TC synapses, the latter of which form of the basis of open loops; outside of this particular functional context, synaptic heterogeneity is broadly speculated to improve the versatility, efficiency, speed, and metabolic economy associated with signal processing (Lengler et al., 2013). If open-loop TC-TRN-TC architectures are indeed to be found within a larger, synaptically diverse thalamoreticular milieu, characterized by variable synaptic divergence and differing densities of sensory, cortical, and other extrinsic innervation, it is moreover reasonable to expect that the degree to which the mode of propagation mediated by these synaptic motifs prevails would vary across thalamic and cortical regions.
To what extent might the functionality of open loops between the dorsal thalamus and TRN depend on arousal state? While sleep and other depressed states of consciousness (e.g., those pharmacologically induced) are associated with thalamic hyperpolarization and therefore tend to amplify both postinhibitory rebound, which underlay the signal propagation in our model, and, by extension, low-threshold bursting McCormick and Bal, 1997;Weyand et al., 2001;Urbain et al., 2019), both phenomena have also been documented in the thalamic relay neurons of awake animals (Guido and Weyand, 1995;Fanselow et al., 2001;Ortun~o et al., 2014), and then sometimes selectively in response to particular stimuli (Lesica and Stanley, 2004;Wang et al., 2007). Thus, while there is no mechanistic basis on which to assume that propagation through open TC-TRN-TC loops would be restricted to a particular state of wakefulness, the widespread modulatory afferents received by both the thalamus and TRN from brain areas including the prefrontal cortex, basal forebrain, amygdala, and brainstem leave little doubt that any form of intrathalamic or cortical signaling supported by these synaptic architectures would be highly state-dependent (Halassa et al., 2014). Both forthcoming physiological investigation and future modeling studies will be able to evaluate such predictions and help provide a full accounting of the role of the various modes of connectivity between brain regions.