Optogenetic stimulation recruits cortical neurons in a morphology-dependent manner

1 Optogenetic stimulation recruits cortical neurons in a 1 morphology-dependent manner 2 3 Short title: Neuron morphology shapes optogenetic response. 4 5 David Berling¹, Luca Baroni¹, Antoine Chaffiol², Gregory Gauvain², Serge Picaud², Ján Antolík¹ 6 1 Faculty of Mathematics and Physics, Charles University, Prague, Czech Republic. 7 2 Institut de la Vision, Sorbonne Université, Paris, France. 8 9 Corresponding author: Ján Antolík [email protected] 10 Figures: 8 11 Tables: 1 12 Abstract: 207 words 13 Introduction: 543 words 14 Discussion: 1487 words 15 Conflicts of interest: Serge Picaud is a founder and consultant for GenSight Biologics and Pixium Vision. 16 Acknowledgements: We thank Thomas Foutz for technical advice on the multi-compartment neuron 17 model and thank Lawrence Humphreys and Jose-David Celdran for helpful discussions on subcellular 18 targeting of opsins. 19 Funding: This work received funding from the European Union’s Horizon 2020 research and innovation 20 programme under grant agreement No. 861423 (Entrain Vision), the Charles University Primus Research 21 Program 20/MED/006 and from the ERDF-Project Brain dynamics, No. CZ.02.01.01/00/22_008/0004643. 22 23 2 Abstract 24 Single-photon optogenetics enables precise, cell-type-specific modulation of neuronal circuits, making it a 25 crucial tool in neuroscience. Its miniaturization in the form of fully implantable wide-field stimulator arrays 26 enables long-term interrogation of cortical circuits and bares promise for Brain-Machine Interfaces for 27 sensory and motor function restoration. However, achieving selective activation of functional cortical 28 representations poses a challenge, as studies show that targeted optogenetic stimulation results in activity 29 spread beyond one functional domain. While recurrent network mechanisms contribute to activity spread, 30 here we demonstrate with detailed simulations of isolated pyramidal neurons from cat of unknown sex that 31 already neuron morphology causes a complex spread of optogenetic activity at the scale of one cortical 32 column. Since the shape of a neuron impacts its optogenetic response, we find that a single stimulator at 33 the cortical surface recruits a complex spatial distribution of neurons that can be inhomogeneous and vary 34 with stimulation intensity and neuronal morphology across layers. We explore strategies to enhance 35 stimulation precision, finding that optimizing stimulator optics may offer more significant improvements than 36 preferentially somatic expression of the opsin through genetic targeting. Our results indicate that, with the 37 right optical setup, single-photon optogenetics can precisely activate isolated neurons at the scale of 38 functional cortical domains spanning several hundred micrometers. 39 40 Significance Statement 41 Sensory features, such as the position or orientation of a visual stimulus, are mapped onto the surface of 42 cortex as functional domains. Their selective activation, that may enable eliciting complex percepts, is 43 intensively pursued for basic science and clinical applications. However, delivery of light into one functional 44 domain in optogenetically transfected cortex results in complex, widespread neuronal activity, spreading 45 beyond the targeted domain. Our computational study reveals that neuron morphology contributes to this 46 diffuse response in a cortical-layer and intensity-dependent manner. We show that enhancing the stimulator 47 optics is more effective than soma-targeting of the opsin in increasing spatial precision of stimulation. Our 48 simulations provide insights for designing optogenetic stimulation protocols and hardware to achieve 49 selective activation of functional domains. 50 4 Introduction 52 Optogenetic stimulation offers neuroscientists a unique tool for manipulating neuronal activity with 53 precision, providing cell-type specificity and high spatial and temporal resolution (Deisseroth, 2015; Emiliani 54 et al., 2022). Until recently, long-term optogenetic wide-field manipulations in higher mammals were limited, 55 but now become feasible through implantable stimulator arrays using single-photon-based stimulation 56 (Rajalingham et al., 2021; Hee Lee et al., 2022). This advancement opens a path towards a sophisticated 57 interrogation of sensory or motor code in meaningful volumes of cortex (Huang et al., 2014; Nassi et al., 58 2015; Chernov et al., 2018; Rajalingham et al., 2021; Wang et al., 2022), as well as development of neuro59 prosthetic systems for vision and hearing loss remediation (Chernov et al., 2018; Kleinlogel et al., 2020; 60 Rajalingham et al., 2021; Sahel et al., 2021). However, these applications require spatially precise 61 stimulation to selectively engage functional cortical representations at a spatial scale of cortical columns 62 (Roe et al., 2020; Seidemann, 2023) requiring sub-millimeter precision (Ohki et al., 2005), which in 63 experiments proved challenging due to the spatial spread of stimulation-evoked responses (Huang et al., 64 2014; Chernov et al., 2018; Wang et al., 2022). The observed spatial spreading and complex distribution 65 of responses to stimulation may be attributed to propagation and transformation of the stimulation-evoked 66 responses in the network (Li et al., 2019; Luboeinski and Tchumatchenko, 2020; Mahrach et al., 2020; 67 Antolík et al., 2021). Alternatively, studies on neuron excitability suggest that optogenetic activation can 68 recruit neurons via distal parts of their morphology (Foutz et al., 2012; Schoeters et al., 2023), raising the 69 question whether direct activation of distant neurons via their spatially extending morphology additionally 70 drives spread-out and spatially complex responses to optogenetic stimulation? 71 Here, we approached this question using computational methods to reveal how optogenetic responses of 72 a layer-2/3 and a layer-5 pyramidal neuron depend on the lateral displacement of an optical fiber placed on 73 the cortical surface. We validated the model's basic optogenetic response dynamics against data from 74 primate retinal ganglion cells (RGCs) recorded across a broad range of stimulation intensities, inducing 75 threshold, sustained and depolarization block responses. 76 Our simulations demonstrate that since distant parts of neuron morphology contribute to activation, 77 optogenetic stimulation directly drives neurons distributed in a complex manner in the cortical network 78 beneath the stimulator. We further find that non-linearities in the optogenetic mechanism and neuronal 79 5 response function cause a stimulation-intensity-dependent redistribution of the activated neurons. Further, 80 we find that due to neuronal morphology, different spatial patterns of cells will be directly optogenetically 81 activated in layer-2/3 than in layer-5. These findings reveal that neuronal morphology causes a diverse 82 stimulation-intensity and layer-dependent spatial distribution of optogenetic responses, that can explain the 83 experimentally observed variation of stimulation outcomes. 84 The simulations predict that surface optogenetic stimulation significantly activates isolated layer-2/3 85 neurons within 100 µm and layer-5 neurons within 200 µm lateral distance to the stimulator, which suggests 86 that surface stimulation may enable selectively recruiting functional neuronal representations encoded at a 87 submillimeter scale of cortical columns. Yet, higher precision may be necessary, since evoked activity could 88 spread beyond the directly activated cortical volume through synaptic transmission in the network. We 89 therefore investigated whether genetics-based or optical improvements to the optogenetic setup could 90 further enhance precision. We simulated soma-targeted channelrhodopsin expression derived from 91 fluorescence measurements in cortical brain slices, finding that spatial precision was only improved for the 92 layer-5 neuron and high channelrhodopsin expression. On the other hand, increasing the collimation of the 93 optical setup improved spatial precision irrespective of the neuron-type and maximized spatial stimulation 94 precision for the tightest light beam explored regardless of whether channelrhodopsin targeting to cell 95 bodies was involved. These results imply that optical improvements of the stimulator are more effective in 96 increasing precision of surface-based stimulation than genetics-based targeting of the opsin to the cell 97 body. 98 99 100 6 Materials and Methods 101 Animals 102 All experiments were conducted in accordance with the National Institutes of Health Guide for the Care and 103 Use of Laboratory Animals. The protocol was approved by the local animal ethics committees and was 104 conducted in accordance with Directive 2010/63/EU of the European Parliament. Cynomolgus macaques 105 (Macaca fascicularis) of foreign origin and unknown sex were used for the study. 106 Experimental Design 107 Experimental data of cell-attached recordings of primate retinal ganglion cells have been partially published 108 previously with a detailed description of methodology (Chaffiol et al., 2017). Electrophysiological data are 109 available from the corresponding author upon reasonable request. 110 AAV production and eye injection 111 Adeno-associated-virus (AAV) backbone plasmid contained human-codon-optimized calcium112 translocating-channelrhodopsin (CatCh) sequence. Gamma-synuclein (SNCG) promoter was cloned into 113 the AAV backbone plasmid. The constructs further contained woodchuck hepatitis virus post-transcriptional 114 regulatory element (WPRE) and bovine growth hormone poly (A). 115 Animals were anesthetized with ketamine/xylazine (10:1 mg/kg). The vitreous was injected with 100 µL of 116 viral vector holding 10¹² viral particles. After injection, the cornea was treated with an ophthalmic steroid 117 and antibiotic ointment. 118 Electrophysiology 119 AAV-treated macaque retinas were recorded in oxygenized (95% O2, 5% CO2) Ames medium (Sigma120 Aldrich), to which a selective group III metabotropic glutamate receptor agonist, L-(+)-2-amino-4121 phosphonobutyric acid (L-AP4, 50 mM, Tocris Bioscience, Bristol, UK) was added. 122 Cell-attached recordings were conducted with an Axon Multiclamp 700B amplifier, using borosilicate glass 123 capillaries (BF100-50-10; Sutter Instruments), pulled to 7-9 MΩ. Retinal ganglion cells were recorded in 124 7 current-clamp configuration (current zero) with electrodes filled with Ames' solution. Retinal ganglion cells 125 were dark-adapted for one hour prior to the recording. 126 Full-field photostimulation was performed with a Polychrome V monochromator (Olympus, Hamburg, 127 Germany), and output light intensities and wavelengths were calibrated. 128 Simulations and computational model 129 Simulation software 130 We simulated neuronal dynamics with the NEURON simulator, version 8.0.0 (Carnevale and Hines, 2006), 131 through its python-interface (Hines et al., 2009) using python 3.6. We used the snakemake workflow 132 manager (Mölder et al., 2021) to organize the simulation workflow, to parallelize, and to manage runs with 133 various parameter sets. 134 Neuron models 135 We simulated neurons with active dendrites using existing multi-compartment models based on 136 reconstructions of a layer-5 and layer-2/3 pyramidal cell (Mainen and Sejnowski, 1996) and 137 electrophysiological dynamics previously identified based on experimental recordings for the layer-5 cell 138 (Hu et al., 2009). Somato-dendritic electrophysiological dynamics included voltage-gated sodium (na) and 139 potassium channels (kv), slow non-inactivating (km), and calcium-dependent potassium channels (kca), as 140 well as high-voltage activated calcium channels (ca). Axonal electrophysiology included lowand high141 threshold sodium channels (na12/na16) as well as voltage-gated potassium channels (kv). Channel 142 conductance throughout the neuron's morphology and reversal potentials are listed in Table 1. Channels 143 are implemented in the corresponding NMODL files (channel abbreviation + suffix ‘.mod’). To model the 144 conductance through Channelrhodopsin-2 channels, we used a NMODL implementation of the channel by 145 Foutz et al. (2012), see “Channelrhodopsin model” for details. As a default, we assumed a uniform density 146 of 130 channels/µm2 throughout the neuron morphology (Nagel et al., 1995; Foutz et al, 2012). This yielded 147 10354945 ChR2 channels for the layer-5 and 2588736 for the layer-2/3 cell (a fourth of the channels in the 148 layer-5 cell). To investigate the effects of preferential expression of ChR2 at the soma (“soma targeting”), 149 we set the somatic density to the otherwise uniformly applied density and decreased ChR density along the 150 neuronal processes according to distributions derived from experimental data (cf. “Derivation of spatial 151 8 channelrhodopsin expression distributions”). In addition to our default expression density of 130 152 channels/µm², we explored a ten-fold stronger expression (1300 channels/µm²) along with the simulations 153 comparing somatic-targeting and no-targeting conditions (Figure 7, 8). 154 Channelrhodopsin model 155 The channelrhodopsin dynamics is based on a four-state Markov model describing the changes in 156 conformation of the molecule (Nagel et al., 2003; Hegemann et al., 2005; Nikolic et al., 2009) and has been 157 previously implemented in NMODL (NEURON) (Foutz et al., 2012). Model parameters are summarized in 158 Foutz et al. (2012) and the original model publication (Nikolic et al., 2009). The model accounts for two 159 open, i.e., ion-conducting, and two closed states. Under light illumination channelrhodopsin can transfer 160 from the closed-1 to the open-1 state. From there it can transfer back to closed-1 or transfer to another 161 open state, open-2, in which its conductance is lower than in open-1. From the open-2 state the molecule 162 can transfer back to open-1 or to closed-2. The closed-2 state can thermally decay into closed-1 or under 163 light illumination photoexcite back into open-2. 164 Optical fiber light model 165 We simulated the propagation of light emitted by an optical fiber inside the cortical tissue accounting for 166 absorption and scattering of photons (absorption coefficient 0.125 mm⁻¹, scattering coefficient 7.37 mm⁻¹ 167 (Foutz et al., 2012; Gradinaru et al., 2009); index of refraction n=1.36 (Vo-Dinh, 2014)) according to the 168 Kubelka-Munk general theory of light propagation (Kubelka and Munk, 1931; Aravanix et al., 2007; Foutz 169 et al., 2012; Vo-Dinh, 2014). The model enables the simulation of various fiber models by adapting the 170 parameters diameter and numerical aperture. The diameter defines the size of the output surface, and the 171 numerical aperture defines the divergence of the emitted light beam. In our simulations we assumed that 172 the optical fiber is placed on top of the cortical surface with its output surface directly connected to the 173 cortical medium. We did not account for light losses involved in the coupling of the light into the cortical 174 tissue. We used a fiber with a 200 µm diameter and numerical aperture of 0.22 as default and explored 175 variations of the diameter of 25-400 µm and of the numerical aperture of 0.1-0.5 specifically [Fig. 7]. 176 Derivation of spatial channelrhodopsin expression distributions 177 9 We derived the spatial distribution of ChR inside the neuronal morphology from experimental data obtained 178 by Shemesh et al. (2017). In this study, mouse cortical neurons were transfected with CoChR (Chloromonas 179 oogama channelrhodopsin) with and without using somatic targeting motifs to restrict the channelrhodopsin 180 expression to the soma. They attached green fluorescent protein (GFP) to CoChR to derive the relative 181 surface density of CoChR through fluorescence measurements, which were taken in 10 µm steps along 182 neurites and normalized to the fluorescence level at the soma. 183 In the condition without somatic targeting, the fluorescence measurements along the neurites had a slight 184 decreasing trend with distance from soma. Considering their uncertainty, they were compatible with a 185 uniform expression throughout the morphology, which we chose as distribution to model the “no targeting” 186 condition [Fig 6A, left]. A uniform distribution in the “no targeting” case is also supported by similar 187 fluorescence measurements for a ChR2-YFP compound in cultured neurons from rats (Grossman et al., 188 2010). 189 Somatic targeting of CoChR resulted in an exponential drop of the relative fluorescence with distance to 190 the soma. However, it was not clear whether the exponential decay plateaued or further decayed beyond 191 their most distant data point at 100 µm distance from soma. To account for these two alternatives, we fitted 192 two distribution models to their soma-targeting data [Fig. 6A, right]: Distribution model 1 represented a 193 relative somatic targeting of ChR. The density decayed exponentially with the distance from soma until it 194 plateaued at about 50 µm distance to soma at a level of about 20% of the somatic density. Since this model 195 represented only an overexpression of ChR at the soma compared to distant dendrites, we referred to it as 196 “soft soma targeting”. Distribution model 2 accounted for a complete decay of the ChR density. It described 197 a similar exponential decay as distribution model 1 but had an additional linear decay which yielded zero 198 ChR density at about 200 µm distance to the soma. We referred to this case as “strict soma targeting”. 199 Analysis of simulated optogenetic responses 200 To obtain spatial profiles of optogenetic responses exhibiting their dependence on the stimulator placement 201 on top of the cortical surface, we simulated optogenetic stimulation with a constant intensity of 200 ms 202 duration. We recorded the membrane voltage at the neuron's soma, interpolated the data at 0.1 ms, and 203 thresholded at 0 mV to count spikes. We converted the spike count to firing rate by dividing by the 204 stimulation duration. We repeated this procedure at stimulator locations starting at the central location 205 10 above the neuron's soma and displacing the stimulator in lateral steps of 25 µm up to 975 µm. We 206 additionally varied the stimulator location along angular direction in 16 steps along the full circle. 207 We decided to compare various conditions requiring the same target response level. We determined the 208 response level as the maximum of the optogenetic response across the stimulator's radial distance 209 coordinate and averaged over its angular coordinate values. Since the response level is dependent on the 210 stimulation intensity, finding a specifically targeted response level required a search algorithm. To find this 211 target response level, we first mapped the optogenetic response for the layer-2/3 and layer-5 neurons on a 212 logarithmic and coarse scale to determine the stimulation intensity thresholds for response and for 213 depolarization block. Using these thresholds as limits, we used a bisection algorithm to find the stimulation 214 intensity resulting in the desired target response level. The algorithm used as a first test value the center 215 between the limits on a logarithmic scale and repeated this procedure while updating the upper/lower limit 216 with the previously tested value depending on whether the resulting response level exceeded or fell below 217 the targeted response level. During this procedure, the target response level was determined from the 218 stimulator-distance-dependent optogenetic response profile only up to a radial distance of 200 µm to 219 accelerate the search. We terminated the optimization procedure when the simulated response level was 220 within 10% tolerance of the targeted one or terminated and discarded the data point if the target response 221 level was not found within 31 steps. 222 Finally, we measured spatial precision of the stimulation by determining the radial stimulator distance at 223 which the optogenetic response across radial stimulator distance and averaged over the stimulator's 224 angular coordinate reached half of the response level. We used linear interpolation if the stimulator 225 distance yielding half of the response level was not covered by sampling. 226 Code Accessibility 227 All original simulation code has been deposited at https://github.com/CSNG-MFF/osmorph/ and is publicly 228 available as of the date of publication. Data and any additional information required to reanalyze the data 229 reported in this paper are available from the corresponding authors upon reasonable request. 230 231 18 cortical surface, we replicated their measured soma-targeted opsin expression distribution in addition to the 394 default uniform distribution representing no targeting in our simulations [Fig. 7A,B; see Methods]. We 395 simulated two alternative somatic targeting conditions as their data were not clear on the presence of 396 residual ChR transfection in the distal dendrites. First, we fitted an exponential decay plus y-offset 397 plateauing at about 20% of the density expressed at the soma (“soft soma targeting”), and second, an 398 exponential and linear decay plus y-offset which decays to zero density at around 200 µm distance from 399 soma (“strict soma targeting”). In addition, we explored the impact of different absolute ChR expression 400 levels in addition to different relative ChR distributions. The default ChR expression density used in our 401 simulations was 130 channels/µm², which matched an empirical estimate of bacteriorhodopsin density in 402 Xenopus oocytes (Nagel et al., 1995). We additionally tested a 10-fold higher expression density (1300 403 channels/µm²). 404 We compared spatial precision between conditions by comparing the response space constant at the same 405 target peak response level. Somatic targeting improved spatial precision only for the layer-5 neuron, and 406 required strict somatic targeting of ChR, as shown for an example target response level of 50 Hz in Figure 407 7C. Higher ChR expression (1300 channels/µm²) further enhanced the precision. As notable precision 408 improvement was restricted to the layer-5 neuron with strict somatic targeting of ChR at high expression 409 density (1300 channels/µm²), we included detailed simulation results only for this condition in Figure 7D. 410 Across conditions (neurons, targeting, ChR expression level), somatic targeting required an increase of 411 stimulation intensity, which amounted to two orders of magnitude in the layer-5 neuron and one order of 412 magnitude in the layer-2/3 neuron [Fig. 7D]. Symmetry of the stimulator-location-dependent response 413 profiles was improved with strict somatic targeting in both neurons [Fig. 7D]. In conclusion, our results 414 suggest that somatic targeting of ChR can improve spatial precision only in the layer-5 neuron, but 415 symmetry in both neurons, at the cost of a strong increase of light irradiance. Its effect further requires zero 416 residual ChR expression in distal dendrites. 417 To understand whether the reported effects of ChR targeting generalized across optical fiber 418 parameterizations, we repeated the simulations for the “no targeting” and “strict targeting” conditions and 419 evaluated the response space constant at a peak response of 50 Hz. Indeed, the substantial improvement 420 of spatial precision through strict targeting of ChR in the layer-5 neuron and the lack of improvement in the 421 19 layer-2/3 neuron were present for most simulated optical fiber parameterizations (fiber diameters d = [25, 422 400] µm, numerical apertures NA = [0.1, 0.5]) [Fig. 8]. 423 Narrowing the stimulator light profile improves precision independently of 424 neuron-type 425 Narrowing the stimulator's light profile decreased the spatial extent of the optogenetic response in both 426 neuron types [Fig. 8]. The relative improvement of their response space constants between the widest 427 (d=400 µm, NA=0.5) and narrowest (d=25 µm, NA=0.1) optical fiber conditions was approximately four-fold 428 for the layer-2/3 neuron and threeor four-fold for the layer-5 neuron, depending on whether somatic 429 targeting was involved or not. With the optical fiber narrowed to the tightest light beam profile explored 430 (d=25 µm, NA=0.1), stimulation precision reached its maximum regardless of whether somatic targeting of 431 channelrhodopsin was involved or not. 432 433 20 Discussion 434 Cortical implants using optogenetic stimulation to engage functionally meaningful neuronal populations are 435 posed to play a key role in unraveling of how cortical circuits process information (Huang et al., 2014; Nassi 436 et al., 2015; Roe et al., 2020; Wang et al., 2022; Seidemann, 2023), develop functional architecture 437 (Mulholland et al., 2024a, 2024b), and in clinical applications as neuro-prosthetic devices for restoration of 438 sensory function, such as vision (Kleinlogel et al., 2020; Antolík et al., 2021; Sahel et al., 2021). For success 439 of many applications, the stimulation must be precise at a scale of hundreds of micrometers, at which many 440 functional neuronal representations, for example orientation preference in primate visual cortex (Ohki et al., 441 2005), are encoded. Our simulations show that optogenetic stimulation directly drives neurons across 442 several hundred micrometers, depending on the shape, size, and cortical depth of neuronal morphology. 443 We find direct optogenetic impact on layer-5 neurons within 250 µm and on layer-2/3 neurons within 100 444 µm distance from the stimulator. These predictions align with experiments in the primary visual cortex of 445 tree shrew, which found significant direct optogenetic activation within 150 µm and none beyond 300 µm 446 when network inputs were abolished via excitatory synaptic blockers (Huang et al., 2014). Crucially, our 447 simulations predict that the impact of neural morphology on direct stimulation precision depends on 448 stimulation intensity in a non-linear fashion, as increased intensity results in both elevated or decreased 449 response depending on secondary factors such as the absolute intensity level and neuron-type, which has 450 practical implications for optogenetic experiments. Finally, while subcellular targeting of channelrhodopsin 451 improves precision under restricted conditions, optical changes to the stimulator enhance precision more 452 generally in our simulations. 453 Our findings can help interpret past and future experiments and guide the design of future optogenetics454 based implants. For example, we find that high light intensity can cause suppression by depolarization block 455 centrally below the stimulator while surrounding areas remain activated [Fig. 3]. This unwanted activation 456 pattern can lead to misinterpretation of experiments and may be difficult to detect. Predicting the intensity 457 triggering depolarization block is challenging due to variable transfection levels between cells and 458 experiments. Hence, when increasing stimulation intensity to recruit more neurons, it is important to monitor 459 neurons in the center to detect suppression by depolarization block. Alternatively, pulsed instead of 460 continuous stimulation could be used to prevent depolarization block (Herman et al., 2014), which we 461 21 observed in our experiments and simulations within 10 ms after stimulation onset, hence requiring pulse 462 frequencies larger 100 Hz to avoid its onset. Interestingly, depolarization block could also be used 463 intentionally to suppress neurons beneath the stimulator while activating those around it. 464 The second implication concerns the spatial distribution of optogenetic excitation, which varies with 465 stimulation intensity in pyramidal neurons of layers 2/3 and 5 [Fig. 4A]. Modulating intensity not only affects 466 response strength but also allows optimization of the distribution of excitation [Fig. 4B,C]. However, 467 achieving central activation of both neuron populations beneath the stimulator requires restricting the 468 intensity to a narrow range, limiting flexibility in controlling response levels. 469 Third, we found that increasing stimulation intensity not only raises the response level but also extends the 470 distance at which a neuron is stimulated. Similarly, optogenetic excitation of inhibitory neurons in mice 471 showed a concomitant increase of the inhibitory neurons' response strength and spatial response extent 472 with higher intensity (Li et al., 2019). Since synaptic activity was not blocked and the observed response 473 extent (> 1 mm) exceeded our findings for isolated neurons, it suggests that network mechanisms may 474 amplify this effect. Therefore, while increasing stimulation intensity amplifies neuronal responses, it also 475 broadens their spatial reach. Vice versa, lowering intensity could lead to more spatially specific responses. 476 This poses important constraints for designing stimulation devices and protocols, which must balance 477 activation strength with precision, for which our simulations could provide guidance. 478 Fourth, we observed that somatic targeting of ChR could more than double the stimulation precision in 479 layer-5 but had minimal effect in layer-2/3 neurons. This improvement was only seen with high ChR 480 expression (1300 channels/µm²) and no residual expression in the apical tuft, highlighting that precision 481 enhancement requires a potent, well-targeted optogenetic construct. However, somatic targeting required 482 up to two orders of magnitude higher light intensities to reach the same response [Fig. 7B,C], raising 483 concerns about phototoxicity in superficial layers. Similarly, narrowing the stimulator's light emission profile 484 to improve precision also raises the risk of phototoxicity as it increases light intensity close to the output 485 [Fig. 8]. Using more light-sensitive opsins may help alleviate this risk by enabling stimulation at lower 486 intensities (Sridharan et al., 2022). 487 Our study predicts direct optogenetic activation at distances greater than 200 µm from the neuron's soma. 488 As shown in Fig. 6, these contributions depend on dendritic excitability and can be suppressed or enhanced. 489 22 While more excitable dendrites do not change the overall pattern of increased response levels and reduced 490 spatial precision with higher stimulation intensity, we find that neurons with less excitable dendrites may 491 offer better precision at intermediate compared to low intensities. Additionally, since lateral responses are 492 suppressed in neurons with less excitable dendrites, this has a similar effect in improving precision to 493 somatic targeting of the opsin. 494 Various computational studies investigated optogenetic stimulation of single neurons (Foutz et al., 2012; 495 Grossman et al., 2011, 2013; Jarvis et al., 2018; Schoeters et al., 2023), but only a few (Foutz et al., 2012; 496 Schoeters et al., 2023) examined how the stimulator's position in relation to the neuron affects the 497 optogenetic response. These studies focused on determining the stimulation intensity required to induce 498 spiking, showing that dendrites contribute to activation and that altering ChR distribution across the cell can 499 lower the stimulation threshold. Our results complement these studies by examining spatial optogenetic 500 responses beyond threshold, demonstrating that response characteristics change with varied stimulation 501 intensity, which is important for the design of spatially coordinated optogenetic interventions. 502 There are several limitations to this study that warrant further research. First, this study focuses on 503 excitatory neurons, providing only limited insight into the response characteristics of inhibitory neurons. Our 504 simulations revealed that proximity to the optogenetic stimulator increases a neuron's sensitivity to 505 stimulation. Given that pyramidal cells which comprise most excitatory neurons typically have an apical 506 dendrite extending toward the cortical surface, these excitatory neurons may, on average, be more sensitive 507 to surface-delivered stimuli than inhibitory neurons located in deeper cortical layers, whose dendrites do 508 not extend in the same way. However, less predictable factors, such as the illuminated ChR-transfected 509 membrane area, intensity threshold for effective ChR-conductance, and specific ChR-transfection levels in 510 inhibitory neurons, can significantly influence their response. Detailed computational studies are needed to 511 explore these variables. 512 Second, since our work only examines single-neuron responses to optogenetic stimulation, it is limited to 513 predictions on the direct optogenetic input to a neuronal network. Network-level effects, like those seen in 514 previous studies showing broadening of activation from random connectivity (Luboeinski and 515 Tchumatchenko, 2020) or functional sharpening in the orientation domain of primary visual cortex (Antolík 516 et al., 2021), must be explored using detailed spiking network simulations. This study paves way towards 517 23 such studies, which can integrate our single-neuron findings with network modeling to account for the 518 impact of neural morphology and optogenetic stimulation at the network level. 519 Third, our model predicts depolarization block at high stimulation intensities, consistent with experimental 520 data, but only in about half of the ganglion cells from a single Macaque retina. The absence of depolarization 521 block in the remaining cells may be attributed to lower opsin expression levels, which could not be 522 experimentally verified. Computationally, we confirmed that lower opsin expression indeed requires higher 523 light intensities to induce depolarization block. However, more experimental evidence is needed to answer 524 whether depolarization block appears more universally at high stimulation intensities across various neuron 525 types. 526 Fourth, to understand how large neuronal populations respond to stimulation, it is essential to characterize 527 neurons with diverse morphologies and biophysical properties beyond the neurons examined here. For 528 example, discovering different response characteristics between excitatory and inhibitory neurons could 529 enable varied stimulation strategies, leveraging their distinct roles in network activation. Furthermore, future 530 work could apply our analysis to computational models of broader variety of the same neuron-type, helping 531 to quantify the effects of inter-neuron variations and potentially detect errors in our model parameterizations, 532 which, while well-established, rely on a limited set of experimentally challenging measurements. 533 Finally, we omitted synaptic background activity, which could interact with the charge introduced through 534 ChR photoactivation, potentially shifting the activation threshold or influencing response characteristics if 535 the synaptic input is spatially biased. Accounting for the spatio-temporal properties of synaptic inputs 536 accurately would require full network simulations, which is beyond the scope of this study. Furthermore, 537 our study considers temporally steady light stimulation. Pulsed stimulation is often used in practice (see 538 Grossman et al. (2011) for a comparison), thus, studying how parameters of such protocol impact the 539 response characteristics is of great interest. Another limitation is that our neuron models did not include 540 axonal arbors. Previous modeling by Foutz et al. (2012) suggests that large axon bundles, such as those 541 in layer 2/3 of the primary visual cortex in cats (Stepanyants et al., 2008), could create additional clusters 542 of stimulator positions that activate the neuron. However, such response clusters have not been observed 543 in optogenetic experiments in the primary visual cortex of tree shrews (Huang et al., 2014). Future work 544 incorporating neuronal reconstructions with axonal arbors could clarify their role in optogenetic activation. 545 25 References 547 Antolik J, Sabatier Q, Galle C, Frégnac Y, Benosman R (2021) Assessment of optogenetically-driven 548 strategies for prosthetic restoration of cortical vision in large-scale neural simulation of V1. Sci Rep 549 11:10783. doi: 10.1038/s41598-021-88960-8 550 Aravanis AM, Wang L-P, Zhang F, Meltzer LA, Mogri MZ, Schneider MB, Deisseroth K (2007) An optical 551 neural interface: in vivo control of rodent motor cortex with integrated fiberoptic and optogenetic technology. 552 J Neural Eng 4:S143. doi: 10.1088/1741-2560/4/3/S02 553 Binzegger T, Douglas RJ, Martin KAC (2004) A Quantitative Map of the Circuit of Cat Primary Visual Cortex. 554 J Neurosci 24:8441–8453. doi: 10.1523/JNEUROSCI.1400-04.2004 555 Carnevale NT, Hines ML (2006) The NEURON Book. Cambridge: Cambridge University Press. 556 Chaffiol A, Caplette R, Jaillard C, Brazhnikova E, Desrosiers M, Dubus E, Duhamel L, Macé E, Marre O, 557 Benoit P, Hantraye P, Bemelmans A-P, Bamberg E, Duebel J, Sahel J-A, Picaud S, Dalkara D (2017) A 558 New Promoter Allows Optogenetic Vision Restoration with Enhanced Sensitivity in Macaque Retina. Mol 559 Ther 25:2546–2560. doi: 10.1016/j.ymthe.2017.07.011 560 Chernov MM, Friedman RM, Chen G, Stoner GR, Roe AW (2018) Functionally specific optogenetic 561 modulation in primate visual cortex. Proc Natl Acad Sci USA 115:10505–10510. doi: 562 10.1073/pnas.1802018115 563 Deisseroth K (2015) Optogenetics: 10 years of microbial opsins in neuroscience. Nat Neurosci 18:1213– 564 1225. doi: 10.1038/nn.4091 565 Emiliani V, Entcheva E, Hedrich R, Hegemann P, Konrad KR, Lüscher C, Mahn M, Pan Z-H, Sims RR, 566 Vierock J, Yizhar O (2022) Optogenetics for light control of biological systems. Nat Rev Methods Primers 567 2:1–25. doi: 10.1038/s43586-022-00136-4 568 Foutz TJ, Arlow RL, McIntyre CC (2012) Theoretical principles underlying optical stimulation of a 569 channelrhodopsin-2 positive pyramidal neuron. J Neurophysiol 107:3235–3245. doi: 570 10.1152/jn.00501.2011 571 26 Gradinaru V, Mogri M, Thompson KR, Henderson JM, Deisseroth K (2009) Optical Deconstruction of 572 Parkinsonian Neural Circuitry. Science 324:354–359. doi: 10.1126/science.1167093 573 Grossman N, Nikolic K, Toumazou C, Degenaar P (2011) Modeling Study of the Light Stimulation of a 574 Neuron Cell With Channelrhodopsin-2 Mutants. IEEE T Bio-Med Eng 58:1742–1751. doi: 575 10.1109/TBME.2011.2114883 576 Grossman N, Poher V, Grubb MS, Kennedy GT, Nikolic K, McGovern B, Palmini RB, Gong Z, Drakakis 577 EM, Neil MAA, Dawson MD, Burrone J, Degenaar P (2010) Multi-site optical excitation using ChR2 and 578 micro-LED array. J Neural Eng 7:016004. doi: 10.1088/1741-2560/7/1/016004 579 Grossman N, Simiaki V, Martinet C, Toumazou C, Schultz SR, Nikolic K (2013) The spatial pattern of light 580 determines the kinetics and modulates backpropagation of optogenetic action potentials. J Comput 581 Neurosci 34:477–488. doi: 10.1007/s10827-012-0431-7 582 Hee Lee J, Lee S, Kim D, Jae Lee K (2022) Implantable Micro-Light-Emitting Diode (µLED)-based 583 optogenetic interfaces toward human applications. Adv Drug Deliver Rev 187:114399. doi: 584 10.1016/j.addr.2022.114399 585 Hegemann P, Ehlenbeck S, Gradmann D (2005) Multiple Photocycles of Channelrhodopsin. Biophys J 586 89:3911–3918. doi: 10.1529/biophysj.105.069716 587 Herman AM, Huang L, Murphey DK, Garcia I, Arenkiel BR (2014) Cell type-specific and time-dependent 588 light exposure contribute to silencing in neurons expressing Channelrhodopsin-2. eLife 3:e01481. doi: 589 10.7554/eLife.01481 590 Hines M, Davison A, Muller E (2009) NEURON and Python. Front Neuroinform 3 591 Hu W, Tian C, Li T, Yang M, Hou H, Shu Y (2009) Distinct contributions of Nav1.6 and Nav1.2 in action 592 potential initiation and backpropagation. Nat Neurosci 12:996–1002. doi: 10.1038/nn.2359 593 Huang X, Elyada YM, Bosking WH, Walker T, Fitzpatrick D (2014) Optogenetic Assessment of Horizontal 594 Interactions in Primary Visual Cortex. J Neurosci 34:4976–4990. doi: 10.1523/JNEUROSCI.4116-13.2014 595 27 Jarvis S, Nikolic K, Schultz SR (2018) Neuronal gain modulability is determined by dendritic morphology: 596 A computational optogenetic study. PLOS Comput Biol 14:e1006027. doi: 10.1371/journal.pcbi.1006027 597 Kleinlogel S, Vogl C, Jeschke M, Neef J, Moser T (2020) Emerging Approaches for Restoration of Hearing 598 and Vision. Physiol Rev 100:1467–1525. doi: 10.1152/physrev.00035.2019 599 Kubelka P, Munk F (1931) Ein beitrag zur optik der farbanstriche. Z Tech Physik 12, 593– 601. 600 Li N, Chen S, Guo ZV, Chen H, Huo Y, Inagaki HK, Chen G, Davis C, Hansel D, Guo C, Svoboda K (2019) 601 Spatiotemporal constraints on optogenetic inactivation in cortical circuits. eLife 8:e48622. doi: 602 10.7554/eLife.48622 603 Luboeinski J, Tchumatchenko T (2020) Nonlinear response characteristics of neural networks and single 604 neurons undergoing optogenetic excitation. Netw Neurosci 4:852–870. doi: 10.1162/netn_a_00154 605 Mahrach A, Chen G, Li N, van Vreeswijk C, Hansel D (2020) Mechanisms underlying the response of 606 mouse cortical networks to optogenetic manipulation. eLife 9:e49967. doi: 10.7554/eLife.49967 607 Mainen ZF, Sejnowski TJ (1996) Influence of dendritic structure on firing pattern in model neocortical 608 neurons. Nature 382:363–366. doi: 10.1038/382363a0 609 Mölder F, Jablonski KP, Letcher B, Hall MB, Tomkins-Tinch CH, Sochat V, Forster J, Lee S, Twardziok SO, 610 Kanitz A, Wilm A, Holtgrewe M, Rahmann S, Nahnsen S, Köster J (2021) Sustainable data analysis with 611 Snakemake. F1000Res 10:33. doi: 10.12688/f1000research.29032.2 612 Mulholland HN, Jayakumar H, Farinella DM, Smith GB (2024a) All-optical interrogation of millimeter-scale 613 networks and application to developing ferret cortex. J Neurosci Meth 403:110051. doi: 614 10.1016/j.jneumeth.2023.110051 615 Mulholland HN, Kaschube M, Smith GB (2024b) Self-organization of modular activity in immature cortical 616 networks. Nat Commun 15:4145. doi: 10.1038/s41467-024-48341-x 617 Nagel G, Möckel B, Büldt G, Bamberg E (1995) Functional expression of bacteriorhodopsin in oocytes 618 allows direct measurement of voltage dependence of light induced H+ pumping. FEBS Lett 377:263–266. 619 doi: 10.1016/0014-5793(95)01356-3 620 28 Nagel G, Szellas T, Huhn W, Kateriya S, Adeishvili N, Berthold P, Ollig D, Hegemann P, Bamberg E (2003) 621 Channelrhodopsin-2, a directly light-gated cation-selective membrane channel. Proc Natl Acad Sci USA 622 100:13940–13945. doi: 10.1073/pnas.1936192100 623 Nassi JJ, Avery MC, Cetin AH, Roe AW, Reynolds JH (2015) Optogenetic Activation of Normalization in 624 Alert Macaque Visual Cortex. Neuron 86:1504–1517. doi: 10.1016/j.neuron.2015.05.040 625 Nikolic K, Grossman N, Grubb MS, Burrone J, Toumazou C, Degenaar P (2009) Photocycles of 626 Channelrhodopsin-2. Photochem Photobiol 85:400–411. doi: 10.1111/j.1751-1097.2008.00460.x 627 Ohki K, Chung S, Ch’ng YH, Kara P, Reid RC (2005) Functional imaging with cellular resolution reveals 628 precise micro-architecture in visual cortex. Nature 433:597–603. doi: 10.1038/nature03274 629 Rajalingham R, Sorenson M, Azadi R, Bohn S, DiCarlo JJ, Afraz A (2021) Chronically implantable LED 630 arrays for behavioral optogenetics in primates. Nat Methods 18:1112–1116. doi: 10.1038/s41592-021631 01238-9 632 Roe AW, Chen G, Xu AG, Hu J (2020) A roadmap to a columnar visual cortical prosthetic. Curr Opin Physiol 633 16:68–78. doi: 10.1016/j.cophys.2020.06.009 634 Sahel J-A, Boulanger-Scemama E, Pagot C, Arleo A, Galluppi F, Martel JN, Esposti SD, Delaux A, de Saint 635 Aubert J-B, de Montleau C, Gutman E, Audo I, Duebel J, Picaud S, Dalkara D, Blouin L, Taiel M, Roska B 636 (2021) Partial recovery of visual function in a blind patient after optogenetic therapy. Nat Med 27:1223– 637 1229. doi: 10.1038/s41591-021-01351-4 638 Schoeters R, Tarnaud T, Weyn L, Joseph W, Raedt R, Tanghe E (2023) Quantitative analysis of the 639 optogenetic excitability of CA1 neurons. Front Comput Neurosci 17, 1229715. 640 https://doi.org/10.3389/fncom.2023.1229715 641 Seidemann E (2023) Invited Session IV: Studies of the visual cortex with sub-millimeter resolution: Toward 642 an all-optical bi-directional interrogation of topographic population codes in primate cortex. J Vis 23:26. doi: 643 10.1167/jov.23.11.26 644 35 712 Figure 5: Comparison of spatial stimulation precision of layer 5 and 2/3 neurons at equal target 713 response levels. (A) Examples of optogenetic response profiles and response space constant at a target 714 response level of 50 Hz (top). Determination of ”peak response“ (maximum) and ”response space constant” 715 (radial fiber offset at half-maximum) from the stimulator-offset dependent optogenetic response. (B) Peak 716 response of layer-2/3 (red) and layer-5 neuron (blue) plotted over stimulation intensity at the stimulator 717 output (top). Response space constant for both neurons plotted over stimulation intensity (bottom). (C) 718 Response space constant plotted against targeted peak response, same color-coding as in (B). The 719 response space constant is not directly dependent on peak response, instead both variables are dependent 720 on the underlying variable stimulation intensity. Diverging of the response space constant at certain peak 721 response level is occurring due to depolarization block at high underlying stimulation intensities. 722 36 723 Figure 6: Impact of dendritic sodium and leak conductance on optogenetic response. (A) 724 Optogenetic responses of layer-5 pyramidal neuron depending on stimulator location for varying dendritic 725 excitability (rows) and different stimulation intensities (columns). Rows represent (1) default, (2) 50% 726 decreased sodium (Na) conductance, (3) 50% increased sodium (Na) conductance, (4) 50% decreased 727 leak conductance, and (5) 50% increased leak conductance. (B) Peak response and (C) response space 728 constant of layer-5-neuron response plotted over stimulation intensity at the stimulator output depending 729 on excitability conditions from (A). (D) Response space constant plotted against targeted peak response, 730 same color-coding as in (B), (C). 731 37 732 Figure 7: Preferential somatic ChR expression improves spatial precision only in the layer-5 733 pyramidal neuron. (A) Normalized channelrhodopsin expression distributions (solid lines) fitted to GFP734 fluorescence measurements (discrete data, mean and std) (34). (B) Expression distributions from (A) 735 applied to the layer-5 pyramidal neuron. (C) Comparison of response space constants observed at 50 Hz 736 target peak response for the layer-5 and layer-2/3 pyramidal neuron type at channelrhodopsin expression 737 levels of 130 channels/µm² and 1300 channels/µm². (D) Detailed results of the effect of somatic targeting 738 for the layer-5 neuron and an expression level of 1300 channels/µm² for which strict somatic targeting 739 achieved the strongest improvements in spatial precision. The 3D plot (upper left sub-panel) illustrates the 740 dependence of peak response and response space constant on stimulation intensity. The bottom projection 741 38 shows the dependence of the peak response on stimulation intensity, which is also plotted in 2D in the 742 upper right sub-panel. Stars in the upper right sub-panel mark the data points at 50 Hz target peak 743 response, from which the bar plot in (C) was generated. Spatial response profiles corresponding to these 744 data points are plotted below (lower right sub-panel). The lower left sub-panel shows the behavior of the 745 response space constant with peak response, which is also plotted as projection in the 3D plot (upper left 746 sub-panel). 747 748 39 749 Figure 8: Narrowing optical fiber diameter and numerical aperture improves stimulation precision. 750 (A) Light intensity levels throughout cortical space with layer-5 and layer-2/3 pyramidal neuron 751 superimposed for several example optical fiber parameterizations. (B) Dependence of the response space 752 constant at a target peak response of 50 Hz on optical fiber diameter and numerical aperture for the layer753 5 and layer-2/3 neuron and using no (red) and strict somatic targeting of ChR2 (blue). 754 755 40 Table 1: Ion channel distribution and properties. Conductance (g) in mS/cm² and reversal potential (E) 756 in mV of the ion channels used in the neuron models (columns) depending on their location in the 757 morphology (rows). 758 na na12 na16 ca kca km kv Soma g 80 - - 0.3 3 0.3 20 E 60 - - 57.7 -90 -90 -90 Dendrites g 60 - - 0.3 3 0.3 10 E 60 - - 49-59 -90 -90 -90 Main axon g - - - - - - 1500 Axon initial segment g - 20-3072 1281920 - - - 2001000 Axon hill g - 2560 - - - - 100 Node of Ranvier g - - 1600 - - - - Myelin g 20 - - - - - - Axon (all) E 60 60 60 - - - -90 759