Quantitative analysis of the optogenetic excitability of CA1 neurons

Introduction: Optogenetics has emerged as a promising technique for modulating neuronal activity and holds potential for the treatment of neurological disorders such as temporal lobe epilepsy (TLE). However, clinical translation still faces many challenges. This in-silico study aims to enhance the understanding of optogenetic excitability in CA1 cells and to identify strategies for improving stimulation protocols.

Methods: Employing state-of-the-art computational models coupled with Monte Carlo simulated light propagation, the optogenetic excitability of four CA1 cells, two pyramidal and two interneurons, expressing ChR2(H134R) is investigated.

Results and discussion: The results demonstrate that confining the opsin to specific neuronal membrane compartments significantly improves excitability. An improvement is also achieved by focusing the light beam on the most excitable cell region. Moreover, the perpendicular orientation of the optical fiber relative to the somato-dendritic axis yields superior results. Inter-cell variability is observed, highlighting the importance of considering neuron degeneracy when designing optogenetic tools. Opsin confinement to the basal dendrites of the pyramidal cells renders the neuron the most excitable. A global sensitivity analysis identified opsin location and expression level as having the greatest impact on simulation outcomes. The error reduction of simulation outcome due to coupling of neuron modeling with light propagation is shown. The results promote spatial confinement and increased opsin expression levels as important improvement strategies. On the other hand, uncertainties in these parameters limit precise determination of the irradiance thresholds. This study provides valuable insights on optogenetic excitability of CA1 cells useful for the development of improved optogenetic stimulation protocols for, for instance, TLE treatment.

1. Introduction

Optogenetics is a technique that can be used to manipulate cellular activity with light (Boyden et al., 2005; Deisseroth, 2011). Specific cell types are sensitized to light via the introduction of gene constructs that encode for optogenetic actuators (Rost et al., 2017). Opsins are light-gated ion channels, pumps or receptors that when illuminated allow hyper- or depolarizing currents through the cell membrane (Deisseroth, 2011). The expression of these opsins is controlled by promoters active in the target cells (Rost et al., 2017). Consequently, when applied to the brain, optogenetics provides the possibility to target specific neuronal populations which makes the technique highly promising to treat a variety of neurological disorders, such as epilepsy (Carrette et al., 2015; White et al., 2020).

Temporal lobe epilepsy (TLE) is one of the most common forms of epilepsy. 30% of patients with TLE cannot be helped with anti-epileptic drugs (Cela and Sjöström, 2019). Based on the premise that temporal lobe seizures originate in the hippocampus, inhibition of the excitatory subpopulations is a possible strategy for TLE seizure suppression (Tønnesen and Kokaia, 2017; Walker and Kullmann, 2020). Using the anion selective opsin halorhodopsin (NpHR) to inhibit hippocampal pyramidal cells has resulted in a decrease of seizure activity (Tønnesen et al., 2009; Krook-Magnuson et al., 2013, 2015). Excitation of the inhibitory hippocampal interneurons using cation conducting opsins, such as channelrhodopsin-2 (ChR2) and variants, has also yielded promising results (Krook-Magnuson et al., 2013; Ladas et al., 2015; Assaf and Schiller, 2016; Lu et al., 2016). Furthermore, optogenetics can be employed in studies on the initiation and propagation of seizures, and to investigate the role different neuronal populations play in the seizure dynamics (Krook-Magnuson et al., 2015; Cela et al., 2019).

These studies are very promising but the development of clinical applications using optogenetics within the brain still faces some challenges (Cela and Sjöström, 2019; White et al., 2020). A first challenge is having a good understanding of the optogenetic effect, because aforementioned seizure-reduction strategies can have ictogenic effects as well (Assaf and Schiller, 2016; Wykes et al., 2016; Tønnesen and Kokaia, 2017). Another challenge is bringing sufficient light into the brain. The optic frequencies have a poor penetration in brain tissue necessitating in-vivo implantation of the illumination source. In order to minimize structural damage due to intrusion, the light source's dimensions are limited (Vandekerckhove et al., 2020). As a result only a small volume is illuminated with a single optical fiber, which could be insufficient for the human brain (Tønnesen and Kokaia, 2017). Furthermore, light absorption by the brain tissue causes heating. Therefore the light intensity should be limited in order to prevent brain damage (Stujenske et al., 2015; Shin and Kwon, 2016; Owen et al., 2019; Peixoto et al., 2020; Acharya et al., 2021).

To minimize stimulation power, one improvement could be the transition to red-shifted opsins (Yizhar et al., 2011; Klapoetke et al., 2014). (Infra)red light gets less absorbed resulting in higher illumination volumes and less heating (Yaroslavsky et al., 2002; Johansson, 2010; Stujenske et al., 2015; Lehtinen et al., 2022). Another solution is increasing the efficiency such that higher photocurrents are generated at lower irradiances. Improving single channel conductance or altering channel kinetics are possible strategies (Berndt et al., 2011; Yizhar et al., 2011; Mattis et al., 2012). Another approach is enhancing membrane expression (Gradinaru et al., 2008, 2010; Rost et al., 2017). A final strategy is spatial confinement of the opsin to specific neuronal membrane compartments (Rost et al., 2017, 2022; Mahn et al., 2018; Messier et al., 2018).

Computational modeling can aid in developing a better understanding of the optogenetic response while avoiding in-vivo animal testing (Williams et al., 2013). In 2006, a first mathematical description of the channel dynamics of ChR2 was published (Nikolic et al., 2006). Since then, several studies have been published that use this and other, more accurate and efficient models in combination with neuron models to design in silico experiments (Nikolic et al., 2009; Grossman et al., 2011; Schneider et al., 2013; Williams et al., 2013; Gupta et al., 2019; Bansal et al., 2021; Schoeters et al., 2021). These experiments allow for easy parameter manipulation and exploration of the stimulation parameter space. Opsin expression can be spatially constrained, expression in different cell types can be tested, and interaction with the 3D light intensity profile can be evaluated. This optic field can be obtained via Monte Carlo based simulations (Wang et al., 1995) that simultaneously can be used to investigate the effect of tissue optical properties (Stujenske et al., 2015; Shin and Kwon, 2016; Acharya et al., 2021).

These optical field studies have shown that the light intensity spatial variation occurs on a neuronal scale. Studies combining optic field simulation with optogenetic response of neurons are scarce. Foutz et al. (2012) investigated this interaction with a multicompartment, neocortical layer V, pyramidal model. Arlow et al. (2013) performed a study on a myelinated axon MRG model. They identified a complex interplay of various simulation parameters on the activation threshold. On the other hand, Grossman et al. (2013) showed that opsin spatial confinement has an impact on action potential generation and propagation. This non-uniform opsin expression has been touched upon by Foutz et al. (2012) as well. In this study, the optogenetic excitability of cornu ammonis (CA1) cells is investigated with the aim of gaining insights that could help guide optogenetic experiments concerning the suppression and initiation of seizures. The novelty is the use of morphologically reconstructed and data-driven biophysical models of CA1 pyramidal neurons and interneurons (Migliore et al., 2018) that are extended to include ChR2(H134R) dynamics, and are subsequently subjected to a Monte Carlo simulated optic field. The effect of opsin expression level and spatial confinement on stimulation thresholds for multiple pulse durations is determined. Furthermore, the impact of various uncertain parameters, e.g., optical field properties, cell to optical fiber orientation and 3D structural cell morphology, is quantified. Based on these results, possible subcellular improvement strategies coupled with ideal optical fiber positioning are identified.

2. Methods and materials

The CA1 cell models used to test the optogenetic excitability are described below. Next, the light intensity fields determined via the Monte Carlo (M.C.) method and the opsin model are elaborated on. Finally, the methods and metrics used to analyze the optogenetic response are explained.

2.1. Neuron models

Four models from Migliore et al. (2018) are used in this study: two models of CA1 pyramidal cells and two CA1 interneurons located in the stratum pyramidale, i.e., a pravalbumin positive basket cell and a bistratified cell. All models have a different three-dimensional structural morphology (see Figure 1A). The cell identification number of the pyramidal cells pyr₁ and pyr₂ is mpg141208_B_idA and mpg141209_A_idA, respectively. The bistratified cell (int₁) has identification number 980513B and the basket cell (int₂) has number 060314AM2. Although the pyramidal cells have a different structural morphology, they are classified under the same morphological type (m-type, Markram et al., 2015; Romani et al., 2022). To avoid confusion, in this manuscript we will talk about difference in cell instead of morphology when comparing pyr₁ vs. pyr₂ (int₁ vs. int₂). Both pyramidal cells are continuous accommodating, while the electrical type or e-type of the interneurons is continuous non accommodating (Petilla convention; Ascoli et al., 2008; Romani et al., 2022). The models include 10 active conductances: a sodium current, three types of potassium (K_DR, K_M, K_A), three types of calcium (Ca_T, Ca_L, Ca_N), two Ca²⁺-dependent potassium currents (K_Ca, K_Cagk), and a nonspecific current. Moreover, a simple calcium accumulation mechanism is included with a single exponential decay of 100 ms. The models were fitted to experimentally obtained voltage traces using a genetic algorithm. For more details of the models and the fitting procedure we refer to Migliore et al. (2018). The models themselves can be found on modelDB (https://modeldb.science/) under accession number 244688. The center of cell's soma is at the origin. The somato-dendritic axis of the neurons is aligned to the z-axis. This somato-dendritic axis is the first principal component determined via a principal component analysis on the 3D-positions of all compartments except the axon. In case of the pyramidal cells, the apical trunk is always directed in the positive z-direction.

FIGURE 1

Figure 1. The simulation framework. (A) Graphical representation of the 3D reconstructed CA1 cell models used in this study: two pyramidal cells (pyr₁ and pyr₂), 1 bistratified (int₁) and 1 basket cell (int₂). The models are adopted from Migliore et al. (2018). This reference frame is applicable throughout the whole study, i.e., the soma is at z = 0 mm, and the somato-dendritic axis is parallel to the z-axis. As depicted, three optrode pitch setups are tested. The color represents the different opsin expression locations. (B) The light intensity profile in gray matter with ${\mu }_a=0.42\text{ m}{\text{m}}^{-1}$, ${\mu }_s=11.33\text{ m}{\text{m}}^{-1}$, and g = 0.88. The irradiance at the cell, with soma 400 μm bellow the fiber, is shown in the inset. (C) The opsin current density at a voltage clamp of −70 mV for increasing irradiance (light to dark: 10^-3 → 10 mW/mm²) of a 1 (left) and 100 ms (right) optical pulse, indicated in gray. (Bottom) The peak (solid) and steady state (dashed) current densities as function of irradiance for a 1 and 100 ms optical pulse in green and blue, respectively. (D) The simulation flowchart with input and uncertain parameters.

2.2. Light field in gray matter

To activate the opsins, the neurons are subjected to light. The light intensity field produced by an optical fiber [100 μm radius, 0.39 numerical aperture (NA); after the optrode in the study of Acharya et al., 2021] is determined via the Monte Carlo method. The used method is based on the direct photon flux recording strategy of Shin and Kwon (2016). The Monte Carlo simulations are solved in cylindrical coordinates for a homogeneous medium as in Stujenske et al. (2015). Because the hippocampus is predominantly a gray matter structure, the in-vivo optical properties of gray matter are selected. The absorption coefficient (μ_a) at 470 nm is obtained by extrapolation of the data reported in Johansson (2010) via a third order polynomial fit on the points between 480 and 550 nm. The scattering coefficient (μ_s) and anisotropy factor (g) are obtained via interpolation from Yaroslavsky et al. (2002). This is often combined into the reduced scattering coefficient $\left[{\mu }_s^{\prime }={\mu }_s(1-g)\right]$. The refractive index (n) is 1.36. The simulations are run with 10⁷ photons and a radial (dr) and axial discretization (dz) of 5 μm. The result of the Monte Carlo simulation with parameters given in Table 1 is shown in Figure 1B, with corresponding irradiance at the neuron level given in the inset. Throughout the study three pitch orientations with respect to the CA1 cells are tested as depicted in Figure 1A.

TABLE 1

Table 1. Parameters used in the Monte Carlo (M.C.) simulations.

2.3. ChR2(H134R) opsin

The selected opsin is a genetically engineered Channelrhodopsin-2 (ChR2) variant: ChR2(H134R). It has an enhanced steady-state photocurrent with a slower closing mechanism than the wildtype ChR2. Its peak operation is at a wavelength of 470 nm (Berndt et al., 2011; Mattis et al., 2012). The opsin's kinetics are modeled with the double two-state opsin model, more specifically the RSRS-final reported in Schoeters et al. (2021), because of its improved computational efficiency compared to alternative four state Markov models. It consists of two independent two-state pairs: O (opening and closing of the channel) and R (change in conductance due to dark-light adaptation). The transmembrane current density (mA/cm²) at a single section is determined as follows:

$\begin{array}{ll}{i}_{\text{ChR}2}=g_{\text{ChR}2}G(V)(O(I_{\text{rr}},V)·R(I_{\text{rr}},V))(V-E_{\text{ChR}2})& (1)\end{array}$

$\begin{array}{ll}G(V)=\frac{1-1.25exp\left(-\frac{V-E_{\text{ChR}2}}{44.52}\right)}{V-E_{\text{ChR}2}}& (2)\end{array}$

$\begin{array}{ll}\frac{dX}{dt}=\frac{X_{\infty }(I_{\text{rr}})-X(t)}{{\tau }_X(I_{\text{rr}},V)}& (3)\end{array}$

$\begin{array}{c}{\tau }_X(I_{\text{rr}},V)={\left[{\tau }_X{(I_{\text{rr}})}^{-1}+{\tau }_X{(V)}^{-1}\right]}^{-1},\text{with }X=O,R\\ {\tau }_O(I_{\text{rr}})=\frac{0.021}{1+exp(1.55)I_{\text{rr}}^{0.37}},\\ {\tau }_O(V)=\frac{23.14}{1+exp(-(V+0.39)/13.19)}\\ {\tau }_R(I_{\text{rr}})=10(1-\frac{0.56}{1+exp(-1.82)I_{\text{rr}}^{-0.50}}-\frac{0.44}{1+exp(17.82)I_{\text{rr}}^{-3.95}})\\ {\tau }_R(V)=\frac{99.74}{1+exp(-(V+38.69)/12.02)}\\ O_{\infty }(I_{\text{rr}})=\frac{1}{1+exp(5.45)I_{\text{rr}}^{-0.70}},\\ R_{\infty }(I_{\text{rr}})=1-\frac{0.77}{1+exp(16.33)I_{\text{rr}}^{-3.62}}\end{array}$

with g_ChR2 the specific conductance (S/cm²), G(V) the rectification function given by Equation 2, V the membrane potential (mV), I_rr the irradiance (mW/mm²), and E_ChR2 = 0 mV the equilibrium potential. The photocurrents for an optical pulse with duration of 1 and 100 ms with increasing irradiance under voltage clamp recording of -70 mV is shown in Figure 1C.

The opsin is spatially confined to specific neuronal membrane compartments (in this manuscript further denoted as subcellular region). In case of the pyramidal cells, it is located in the soma, or distributed over the axon, basal dendrites (basal), apical dendrites (apic), all dendrites (dend = basal∪apic), or all sections (allsec = dend∪soma∪axon). In case of the interneurons, no distinction between apical and basal dendrites is made. The different regions are illustrated via the color code in Figure 1A. In each simulation, g_ChR2 is uniformly spread over the subcellular region of interest. Its value is calculated from a preset total maximum conductance:

$\begin{array}{ll}{G}_{\text{max}}=g_{\text{ChR}2,k}·A_k& (4)\end{array}$

with A_k the total membrane surface area of a subcellular region k. The surfaces are summarized in Table 2. G_max are eight points spaced evenly on a log scale between 10⁻¹ and 10^1.5 μS (for the uniform field nine additional points are included between 10⁻¹ and 10²).

TABLE 2

Table 2. Total membrane surface area of subcellular regions in μm².

Finally, the total temporal averaged current at the excitation threshold (TAC; Williams and Entcheva, 2015) is calculated by:

$\begin{array}{ll}\text{T}\text{A}\text{C}=\frac{1}{\text{pd}}{\displaystyle \sum _{jϵk}}{A}_j{\displaystyle {\int }_{t_0}^{T_{\text{end}}}}{i}_{\text{ChR}2,j}(t)\text{ d}t& (5)\end{array}$

where pd is the pulse duration, t₀ is the pulse onset time, T_end = max(500 ms, t₀+pd+100 ms) (to ensure to include channel closure) and j is every compartmental segment in region k. The pd varied over 5 points logarithmic evenly spaced between 0 and 100 ms (9 pds ϵ [10⁻¹, 10³] are used in the simulations with uniform light field).

2.4. Analyses

The tests and metrics used to analyze the optogenetic response are described below.

2.4.1. Surface of fiber positions for the activation of neurons

As metric to assess the optogenetic excitability, the surface of fiber positions for the activation of neurons (SoFPAN) is determined. This metric is similar to the more well-known volume of tissue activation (VTA; Butson et al., 2007; Duffley et al., 2019). In VTA, the frame of reference is the position of the stimulation device, while here, the frame of reference is the stratum pyramidale with the center of cell's soma at origin. In other words, SoFPAN encompasses the positions at which the optical fiber can be located in order to activate the neuron of interest. This is chosen because of the layered structure of the hippocampus, constraining the cell bodies of the considered neurons to the stratum pyramidale. Because the optical fiber is limited to the plane shown in Figure 1A to reduce the number of simulations, only a 2D plane is explored resulting in a surface instead of a volume being determined.

The intensity threshold (I_th) to elicit an action potential (V>−10 mV) recorded in the soma for a given pulse duration and fiber position is determined first. This threshold value is the intensity at the fiber surface (I_fiber). It is obtained with a titration process using the bisection method for seven iterations (i.e., c_i = (a_i+b_i)/2, with b₀/a₀ = 10). One hundred and twenty-one fiber positions are evaluated with z taking on eleven linearly spaced values in [−400, 700] μm. For a pitch of π/2, x ϵ [−1,000, 4,000] μm with Δx = 500 μm. For the fiber pitch of 0 and π, x ϵ [0, 2,500] μm with Δx = 250 μm. For a given pulse duration, the SoFPAN is then determined as the union of the spatial points for which the threshold is lower than the fiber intensity (true: I_th<I_fiber) multiplied by the discretization surface (i.e., 110 × 250 or 110 × 500 μm² depending on the fiber pitch). A lower bound is determined by counting only the discretization surfaces enclosed by four true-points. In case of the upper bound, the true-false field is dilated first with a 3 × 3 mask before the enclosed-by-four-true-points regions are counted. For the pitches of 0 and π, the obtained SoFPAN is multiplied by two, such that the maximal SoFPAN for all pitches is equal to 5.5 mm². The SoFPAN is calculated for nine logarithmic evenly spaced fiber intensities 0.1 and 1,000 mW/mm².

If the cell's dimensions are assumed to be negligible with respect to the spatial variation of the intensity field, the irradiance is uniform over the whole neuron. The cell can thus be represented as a single point in space. The estimated SoFPAN under this assumption is denoted as SoFPAN_uniform.

2.4.2. Wilcoxon signed-rank test

A single-sided, paired Wilcoxon signed-rank test is used to test whether the results from two populations are significantly different. For the uniform field 20 classes (cell and opsin location: 4× {all, axon, soma, dend} + 2× {basal, apic}) are mutually compared. In case of the M.C. field there are 42 classes (pitch, cell, and opsin location: 3× [4× {all, axon, soma} + 2× {basal}]). A scoring system is used to identify the most excitable setup. Here, a 1, 0.1, and 0.01 is given if the SoFPAN is statistically greater with p-value below 0.001, 0.01, and 0.05, respectively. The maximum score is consequently 41.

2.4.3. Regression

A two step fitting procedure is performed to analyze the effect of pulse duration and expression level on the intensity threshold as suggested by Williams and Entcheva (2015). First Lapicque's formulation is fit to the TAC. This gives the strength-duration relationship like with electrical direct current stimulation (Reilly, 1998). Second a linear regression is performed with independent variables the log₁₀ of $\widehat{\text{TAC}}$ (μA) and G_max (μS).

$\begin{array}{ll}\widehat{\text{TAC}}=\frac{\text{TA}{\text{C}}_0}{1-exp(-\text{pd}/{\tau }_{\text{TAC}})}& (6)\end{array}$

$\begin{array}{ll}{log}_{10}\widehat{I_{\text{th}}}=a_G{log}_{10}{G}_{\text{max}}+a_{\text{pd}}{log}_{10}\widehat{\text{TAC}}+c& (7)\end{array}$

Here, TAC₀, τ_TAC, a_G, a_pd, and c are the regression parameters. ${\hat{I}}_{\text{th}}$ and $\widehat{\text{TAC}}$ are the estimators of the threshold intensity and total averaged current, respectively.

2.4.4. Optimal and worst fiber position

The optimal and worst fiber positions for a given pulse duration are defined as the z-position where the depth of activation (amount of positions along the x-axis with I_th < I_fiber) is the highest or the lowest, respectively. As a tie-breaker, the z-position is selected where the average of TAC along x is, respectively, the lowest or the highest.

2.4.5. Elementary effects

The inputs to the simulations can be divided into two categories (Figure 1D). On the one hand, there are the parameters that are known (e.g., optrode radius, NA, and pitch) or controlled by the user (e.g., pd and I_fiber) in an optogenetic experiment. On the other hand, there are variables that cannot be controlled and contain some uncertainty, e.g., the tissue's optical properties (μ_a, ${\mu }_s^{\prime }$), the neuron's structural morphology (cell: pyr/int₁ vs. pyr/int₂) and orientation (roll, yaw), and the opsin expression [level (G_max) and location]. Because of the expected non-linear interactions and high computational time (±12 h for a 121 position sweep with five pulse durations), the elementary effects method is adopted for global sensitivity analysis (Saltelli et al., 2019).

Six influential factors are investigated, i.e., μ_a, ${\mu }_s^{\prime }$, G_max, cell, opsin location, and roll. The measure used (μ^*) is the mean of the absolute value of r = 16 elementary effects. Also, the standard deviation of the elementary effects (σ) is calculated to track interactions and non-linear effects. The r repetitions are sampled with the radial-based design using Sobol's numbers (Campolongo et al., 2011). μ_a, ${\mu }_s^{\prime }$, and G_max are normally distributed with mean and coefficient of variation (C.V.) given in Table 1. The mean and C.V. of G_max are 1 μS and 0.15. The roll is uniformly distributed between [-π, π[ radians. Cell and location are discretely uniformly distributed, with {pyr₁, pyr₂} ({int₁, int₂}) and {all, axon, soma, basal} ({all, axon, soma}) classes, respectively, for pyramidal (inter-) neurons.

2.5. Software

Simulations were done with NEURON 8.0.0 (Carnevale and Hines, 2006) and Python 3.9.12. on the HPC system with AMD Epyc 7552 processing units, provided by the Flemish Supercomputer Center.

3. Results

The results of this study can be subdivided into three sections. First, there is the optogenetic response under a uniform field. This means that all cell sections receive the same irradiance. This facilitates the analysis of the importance of subcellular regions, pulse durations, and expression levels. Next, the Monte Carlo simulated light field is included. This allows us to indicate the effect of light propagation on optogenetic excitability and provides information on optimal or worst fiber positioning. Finally, an elementary effects study is performed to identify the most sensitive parameters.

3.1. The optogenetic excitability in a uniform light field

Figure 2A depicts the intensity threshold (I_th) and temporal averaged current (TAC) curves as function of pulse duration (pd) of several subcellular and cell combinations. Both decrease with increasing pd until a saturation value is reached. For the same subcellular region (i.e., allsec), the threshold curve is shifted down in case of pyr₁ vs. int₂ (blue solid vs. blue dashed-dotted line). Staying within the same cell, it can be seen that changing the subcellular region (all → basal or soma) can result in a downward shift, as well. The same effects can be seen on the TAC, but less pronounced. Increasing the expression level (G_max) from 1 to 10 μS (blue solid line with circle markers) results in an expected decrease of the threshold. On the other hand, the TAC remains the same. This is highlighted in Figure 2B where the TAC is constant for all pds for a G_max > 1 μS in the allsec-pyr₁ setup. A linear relationship for I_th can be observed for G_max≥0.518 μS. No threshold could be found for G_max < 0.215 μS. This is because the photocurrent saturates at high intensities (see Figure 1C) and therefore cannot compensate for the low G_max. At high pulse durations a decrease in TAC is already observed at higher G_max (< 1 μS). For the corresponding I_th, the photocurrent is biphasic with a peak and steady-state value. At these pulse durations (100 and 1,000 ms) the action potential (AP) occurs during the inactivation from peak to steady state. The occurrence of the peak is much more efficient in eliciting an AP. A large part of the total photocurrent occurs before the AP while this diminishes at higher G_max with monophasic photocurrents. Additionally, during the AP, the current drops to zero due to channel shunting. Consequently, the net effect of this drop on the TAC is larger here than in case of no biphasic photocurrents (at low irradiances) or lower pulse durations, where the AP is elicited during the deactivation phase of the photocurrent.

FIGURE 2

Figure 2. The optogenetic response under a uniform light field. (A) The intensity threshold (I_th) and corresponding total temporal averaged current (TAC, top and bottom, respectively) as function of pulse duration. All lines shown are for a G_max = 1 μS except the circle markers (G_max = 10 μS). (B) The intensity threshold and corresponding total temporal averaged current (top and bottom, respectively) as function of G_max for allsec-pyr₁. (C) The p-value of the mutually, paired Wilcoxon signed-rank tests. Single-sided test if class at the row is lower than class in the column for the I_th (left) and TAC (right). The color code at the top and left indicate the cell model. The cell opsin location combinations are sorted on their excitability score (see Section 2.4.2). (D) The results of the two step regression. The goodness of fit is indicated with the adjusted R-squared measure: ${\bar{R}}_{\text{TAC}}^2$ for the Lapicque fit to the TAC and ${\bar{R}}_{\text{tot}}^2$ for the linear regression fit to the threshold. Box plots are shown for ${\bar{R}}^2$-values and regression parameters, calculated for the different cell-types and opsin locations. (E) Summary of the input impedance at 0 Hz, (left) the local impedance of each section in pyr₁, (right) the summarized impedances of the four tested cells. The gray bars are 100 μm.

The I_th is determined for all pyramidal and continuous non-accomodating interneuron models of Migliore et al. (2018), along with virtual clones optimized using HippoUnit (Saray et al., 2021). The result of these additional 38 pyramidal and 39 interneuron cell models is shown in Supplementary Figure 1. Overall, for a given subcellular opsin location, the variation is within one order of magnitude, except when the opsin is expressed in the axon. This observation holds for both pulse durations of 1 and 100 ms and maximum total conductances of 1 and 10 μS. The selected cells of this manuscript do not appear as outliers except for the axon of the bistratified cell. Therefore, it can be concluded that the selected cells are representative for the population.

Mutually paired Wilcoxon singed-rank tests are performed to identify the most excitable subcellular region. The p-values of a single-sided test are shown in Figure 2C of the threshold intensity and TAC on the left and right, respectively. Each compared population consists of 153 points (9 pds and 17 G_max values). The basal-pyr₁ is the most excitable with a I_th significantly (p < 0.001) lower than any other location-cell combination. The same is true for the TAC. The pyramidal cells are more excitable than the interneurons with pyr₁ the most excitable (p < 0.001). Aside from the axon in pyr₂, the order of the excitability scores is the same in both cells (see Section 2.4.2). The apical dendrites is the least excitable subcellular region of the pyramidal cells (p < 0.001). On average a lower I_th is required for the basket cell than the bistratified cell (p < 0.001). Overall a significantly lower threshold is accompanied with a significantly lower TAC. Some exceptions exist where the TAC is significantly greater (e.g., all-pyr₁ vs. soma-int₂). The median of the relative change in I_th of the location-cell combinations as ranked in Figure 2C is shown in Supplementary Figure 2A. The average of the medians is −22.76%. Therefore, the I_th drops every time on average with 22.76% when moving from the least (axon int₁) to the most (basal pyr₁) excitable location-cell combination. Optimal subcellular expression in a cell can result in I_th drops of >75%. Compared to no subcellular specificity (opsin over the whole cell), I_th reductions of >60% can be achieved in the pyramidal cells, increasing toward 83 and 92% for the basket and bistratified cell, respectively (cf. Table 3). The full analysis of the median relative differences of all different subcellular-cell combinations is shown in Supplementary Figure 3A.

TABLE 3

Table 3. Influence of subcellular opsin expression.

The two step regression is performed on each location-cell combination separately. The combined results are shown in Figure 2D. The variability on TAC is captured by Lapicque's formulation resulting in a median adjusted coefficient of determination (${\overline{R}}_{\text{TAC}}^2$) value of 0.99978. The median rheobase (TAC₀) is 3.17 nA and the median time constant (τ_TAC) is 42.95 ms. The variability of I_th is well-explained by the two step regression model (median ${\overline{R}}_{\text{tot}}^2$ = 0.96388). The median parameter values of G_max and pd (a_G and a_pd) are −1.47 and 0.53, respectively. Even though the latter is positive, I_th decreases with pd due to the Lapicque's formulation. Based on these values, G_max has thus a stronger impact on I_th than the pd.

The input impedance of the cell's sections, measured after 100 ms after initialization at −70 mV, are shown in Figure 2E. The impedance is a proportionality factor on the input current resulting in a given voltage change. Therefore, the voltage change is higher for a higher input impedance given a constant input current. The axons have on average the highest impedance but also the highest spread with many outliers. At the left, it can be seen that the impedance increases from axon hillock to the synapses. The soma has a low impedance compared to the median of the other subcellular regions.

3.2. The optogenetic response in the Monte Carlo light field

In this section, the cells are subjected to an intensity field produced by a 100 μm fiber. The 3D-profile is obtained via the Monte Carlo method in homogeneous gray matter tissue with the default optical parameters as summarized in Table 1. Unlike in the previous section, the irradiance at the different neuron sections will depend now on their relative position with respect to the fiber and its orientation. Consequently, I_th now depends on the position in the reference frame as indicated in Figure 1A. An example is shown in Figure 3A, with the I_th map on the left and corresponding TAC on the right, for an allsec-pyr₁ with G_max = 1.18 μS, pd = 10 ms, and a fiber pitch of π/2. A hotspot is observed around [x = 0.43 mm, z = −0.036 mm], where the threshold is the lowest. Even when the neuron lies behind the fiber, i.e., x < 0 mm, it is still possible to excite the cell but at higher I_th. Two islands can be observed for the TAC with a factor two difference between minimum and maximum. The farther away from the cell the lower the variation in TAC.

FIGURE 3

Figure 3. The optogenetic response in a homogeneous Monte Carlo simulated light field. (A) The threshold intensities (left) of the pyr₁ cell, with opsin distributed over all sections, and corresponding total temporal averaged current (right). In the colormaps, the position of the optical fiber is varied with respect to the soma, which is fixed at (x, z) = (0, 0). (B–D) The Surface of Fiber Positions for the Activation of Neurons, the optimal and worst z-position for neuron activation, respectively, as function of fiber intensity. The shaded area in (B) is enclosed by the upper and lower bound. The lines are the result of the pyr₁ cell with opsin location shown in legend at the right top corner. G_max = 1.179 μS and pd = 10 ms in (A–D). In (C, D) the scatter plots overlay the density, averaged over all cell, location, G_max and pd combinations. (E) Excitability score based on paired Wilcoxon signed-rank tests. Each population consists of 360 (5 pd × 8 G_max × 9 I_fiber values) combinations. The pitch of the fiber is illustrated by the orientation of the fiber icon on the left and valid for the whole row. (F) Relative error of SoFPAN, if the light intensity is considered uniform over the neuron: rel. error = (SoFPAN_uniform—SoFPAN_M.C.)/SoFPAN_M.C., the outliers are not shown. The errors are calculated for all the G_max, pd, pitch, and I_fiber combinations.

Based on these threshold maps, the surface of fiber positions for the activation of neurons (SoFPAN) can be estimated. This is shown in Figure 3B for the three investigated fiber pitches. The SoFPAN is a measure of the excitability of the subcellular optogenetic configuration in the light field. The uncertainty caused by the discretization on the SoFPAN is indicated by the shaded area (see Section 2.4.1). There is a larger uncertainty at lower intensities due to the rough simulation grid. The soma appears to be more excitable than the opsin in all sections (red vs. blue) with already a non-zero SoFPAN for a I_fiber of 0.1 mW/mm². This is true for all fiber pitches. The SoFPAN saturates due to the finite simulation domain.

An optimal and worst z-position for each fiber pitch can be determined, as well. The result for the two cases as above are shown in Figures 3C, D. It can be seen that the positions vary with increasing I_fiber. At low intensities, the optimal position is near the subcellular region (cf. soma, red) or the most excitable region (cf. allsec, blue), which is the basal dendrites. At higher intensities, it is better to illuminate the whole simulation domain, explaining the optimal positions at −0.4 and 0.7 mm for the 0 and π pitches, respectively. Once the whole simulation field is excited, the optimal position is completely determined by the average TAC along the x-axis. This can be seen in the sudden changes at the highest I_fiber (>100mW/mm²). In green the density distribution of all 560 combinations (pd × G_max × loc × cell) is displayed. The optimal position is more concentrated with a 0 or π/2 pitch. On the other hand, the worst position is more concentrated with the π/2 and π pitches. The optimal and worst fiber positions are most distant with the π/2 pitch position.

The excitability of the subcellular region is summarized in Figure 3E. The scores are based on the p-values of the single-sided, paired Wilcoxon signed-rank test, as described in Section 2.4.2. At each pitch individually, the subcellular excitability order remains, generally, the same as under the uniform field (cf. Supplementary Figure 4. Only the axons appear to be susceptible to the fiber pitch position. For instance, the axon-pyr₁ combination has an increased (decreased) excitability under the π (0) fiber pitch. Also at pitch 0, the axon of pyr₂ has a lowered excitability with a SoFPAN significantly lower (p < 0.001) than soma-pyr₁. On the other hand, axon-int₂ has a higher SoFPAN than all-int₂ under the π/2 and π pitches. Basal-pyr₁ is under all pitches the most excitable combination, with a score of 41 at pitch π/2. Overall, pitch π/2 is significantly more excitable than π, which is in turn more excitable than 0 (p < 0.001). The axon-pyr₂, however, has a significantly higher SoFPAN (p < 0.001) under pitch π than π/2. The median of the relative change in SoFPAN of the location-cell-pitch combinations according to this ranking, i.e., lowest to highest score, is shown in Supplementary Figure 2C. The average of medians is now 10.45%. Therefore, the SoFPAN increases every time on average with 10.45% percent when moving from least (axon int₁, pitch = 0) to most (basal pyr₁, pitch = π/2) excitable location-cell-pitch combination. The same analysis is performed when restricted to a single pitch. The average of medians of the relative changes in SoFPAN increase toward 56.86, 25.89, 59.18% for a pitch = 0, π/2, and π, respectively. Matching the ideal fiber location to the subcellular expression of a single cell can result in doubling of the SoFPAN compared to the worst combination. Even a tenfold increase is observed in case of the bistratified cell. Compared to subcellular unspecificity (all), SoFPAN increases of 33–100% can be obtained by specifying subcellular expression (cf. Table 3). The full analysis of the median relative differences of all different subcellular-pitch combinations for each cell separately is shown in Supplementary Figure 3B.

The SoFPAN for a light field that is considered to be uniform over the whole cell (SoFPAN_uniform, see bottom Section 2.4.1) is calculated, as well. The relative error compared to the SoFPAN under a Monte Carlo field, i.e., (SoFPAN_uniform-SoFPAN_M.C.)/SoFPAN_M.C., is shown in Figure 3F. The errors is calculated for all the G_max, pd, pitch, and I_fiber combinations. The error of the soma is negligible. The outliers are not shown. Of these, only 0.85% produce a relative error above 5%. These results validate the method as the soma should not depend on the M.C. field as its section is only one point in the 3D-space. The SoFPAN of basal dendrites gets overestimated, with 18.17% having a relative error above 20%. In case of the pyramidal cells and with the opsin distributed over the whole cell, the SoFPAN is predominantly underestimated and the interquartile range (IQR) of the relative error is 60% of the IQR for the basal dendrites. The estimation for the axon of the int₁ cell is the worst with a median relative error of −50%. At least in one test case of each subcellular region the error is either −100% or infinity (see Table 4).

TABLE 4

Table 4. Summary relative error SoFPAN Monte Carlo vs. SoFPAN uniform (all cells pooled together).

3.3. Parameter uncertainties

There is uncertainty on various parameters used in this study. The optical parameters depend on multiple factors, e.g., tissue and wavelength. The absorption coefficient is extrapolated which introduces an uncertainty, as well. Moreover, the values of gray matter are used while different gradations exist. The effect of a change in optical parameters on the light field in homogeneous tissue, is illustrated in Figure 4A. On the left, the field as used in the section above is depicted. The effect of an increased absorption and decreased reduced scattering coefficient is shown in the middle. The result is a more conical field with higher degradation. A more round field is obtained when the reduced scattering coefficient is higher (cf. right). Also at the cell level, there are multiple sources of uncertainty. In experimental setting, the opsin expression G_max will not be exactly known. Moreover, the subcellular location will probably not be discrete as used in these simulations. Finally, the morphology of the tested cells and its orientation (cf. roll) are fixed. To address the impact of these uncertainties on the output, a global sensitivity analysis is performed. The used approach is a screening method: the Elementary Effects test.

FIGURE 4

Figure 4. Results of the elementary effects study. (A) Intensity fields with different optical parameters; from left to right ${\mu }_a=0.42,0.52,0.35\text{ m}{\text{m}}^{-1}$ and ${\mu }_s^{\prime }=1.36,0.93,2.06$ mm⁻¹. (B, D) The influential parameters on SoFPAN and optimal fiber position, respectively, ranked for 5 pds and 9 fiber intensities (in %). The fiber pitch and cell type are shown on top and left, respectively. (C) The two measures of the elementary effects of three pyramidal cell setups, i.e., pitch = π, π/2, and 0, and I_fiber = 1, 1,000, and 1,000 mW/mm² indicated by circle, diamond and cross, respectively, for a pd = 10 ms. The legend between (B, D) also applies on (C). (E) Normalized optimal z-position over all results with allsec subcellular location of the elementary effects study.

The influence of these six parameters, i.e., cell roll, ${\mu }_s^{\prime }$, μ_a, opsin subcellular location (loc), cell model (cell), and G_max, on the SoFPAN for the three fiber pitches and two cell types (pyr and int) is investigated. For each fiber pitch and cell type, the elementary effects test (EET) is repeated for 5 pd and 9 I_fiber combinations (these 45 combinations correspond with 100% in Figures 4B, D). The rank according to the μ^*-measure is summarized in Figure 4B. For certain pd and I_fiber no differentiation could be made based on the μ^*-measure, explaining the bar height >100% at rank 0. The subcellular location has most frequently the highest impact on the excitability for all six test cases. This is followed by G_max in second place (rank 4). The roll and μ_a have on average the lowest impact. However, in case of the pyramidal cells and π fiber pitch, the roll has occupied the highest rank for some (pd, I_fiber) combinations. On the other hand, ${\mu }_s^{\prime }$ is more important when the fiber pitch is 0 for the pyramidal cells. The μ^* and σ measures of three cases of the pyramidal cells with pd of 10 ms are shown in Figure 4C. The circles indicate the setup where the roll has the highest rank, i.e., pyramidal cell with π pitch and 1 mW/mm². It can be seen that even though it has the highest rank, its measures are in the same range as the other two cases. The diamonds represent the EET of a π/2 pitch at I_fiber = 1, 000 mW/mm². Here, the reduced scattering coefficient has the highest impact. The non-linear and interaction effects are higher in this case, reflected by the higher σ. For the 0-pitch with I_fiber = 1, 000 mW/mm² (cross), the location is ranked highest with a clear difference in μ^*. The effect of cell appears to be more linear than the other parameters indicated by its shift toward lower σ values.

The same analysis is performed on the optimal position. The rank is summarized in Figure 4D. For the interneurons, the subcellular location stays dominant. There is, however, a clear shift toward the optical parameters for the pyramidal cells. Here, the reduced scattering and absorption coefficient are most often ranked highest for the π and 0 pitch, respectively. Also, in case of the π/2 pitch, the absorption coefficient is more important for the optimal position than it is for the SoFPAN. The cell and roll appear to have the lowest effect. While, on the other hand, the cell is important in case of the interneurons. The optimal position for both cell types normalized over all simulations with allsec subcellular location of the EET study, is shown in Figure 4E. These exclude the preset subcellular selectivity. For the interneurons this is more smeared out and a focus toward the stratum pyramidale is observed. On the other hand, for the pyramidal cells there is a clear preference for a position such that the majority of the light reaches the stratum oriens region. At pitch π this is more smeared out due to the possibility to retract (z more positive) the fiber at higher intensities to illuminate a bigger region in the stratum oriens. This is limited for the 0 pitch explaining the high peak at −0.4 mm.

4. Discussion

In this study we focused on the optogenetic excitability of CA1 cells. We attempted to not only gain more insight into the effect of the various stimulation and uncertain parameters but also to identify strategies for increased optogenetic efficiency. These insights are of interest for the development of better stimulation protocols that can be used as treatment for TLE. A broad view is adopted where both the excitability of pyramidal and interneurons is investigated. Even though excitation of inhibitory neurons is one of the two main investigated strategies as treatment of TLE, insights in the excitability of pyramidal neurons can be of interest as well. Like with electrical deep brain stimulation, the latter could be used as counter-irritation (Carrette et al., 2015) with various modes of action (Sprengers et al., 2020) that can be tested. Moreover, stimulation of both types could be beneficial for restoring the excitation-inhibition balance (Wykes et al., 2016).

4.1. Excitability of spatially confined opsin expression

The results show that the optogenetic excitability of CA1 cells depends on various parameters. The irradiance threshold ranges over multiple order of magnitudes. As expected, an increase in expression level (G_max) or pulse duration (pd) results in a decrease of the intensity threshold (I_th). There is also a clear dependence on the subcellular region of opsin expression and variance among different cells. There is no single explanation for the relative excitability of the considered subcellular opsin locations, due to the complex interplay of many non-linear relationships. By comparing the membrane areas and impedances (cf. Table 2 and Figure 2E) of the subcellular regions some observations can be made. For a fixed G_max, the specific channel conductance (g_ChR2) is locally higher for regions with a lower total surface area. Thus, for the same I_fiber, there will locally be a higher depolarizing current to elicit an action potential (AP). This combined with the fact that the AP is measured in the soma (therefore does not have to travel through the cell), explains why the soma-confinement is highly ranked in each cell. The observation concerning the locally raised channel conductance also holds when comparing basal dendrites with apical, all dendrites, and all sections. An argument for why confinement to the basal dendrites is more excitable than the soma could be found by comparing their input impedances. For a fixed depolarizing current, a higher impedance results in a larger membrane depolarization. Because the impedance of the basal dendrites is significantly higher than that of the soma (log-scale in Figure 2E) it will facilitate AP initiation. However, this contradicts the rank of axon-pyr₁. Finally, there is the channel distribution inside the cell itself. The ratio of depolarizing (e.g., Na⁺ and Ca²⁺) and hyperpolarizing (e.g., K⁺) channels defines the membrane threshold for AP initiation. This ratio is double in the axon of pyr₂ compared to pyr₁, motivating the low rank of axon-pyr₁. These observations are in agreement with the findings of Foutz et al. (2012).

Confinement to the basal dendrites of pyr₁ is the most excitable, while the highest I_th is required for the axon of int₁ (cf. Figure 2C). Similar ranking is observed in the SoFPAN calculations (cf. Figure 3E). To identify the effect of endogenous channel distribution, the channel distributions of pyr₁ were imposed on pyr₂ (pyr₃) and vice versa (pyr₄). The mean relative difference of all pd and G_max combinations is shown in Table 5. The I_th of pyr₃ drops with ~30 and ~50% vs. pyr₁ and pyr₂, respectively, for all subcellular regions except for the axon (both cases) and for the basal dendrites in pyr₃ vs. pyr₁. Contrarily, the excitability of pyr₄ drops (higher I_th), except for the basal dendrites for pyr₄ vs pyr₂ or axon in pyr₄ vs. pyr₁. Switching 3D structural morphology while maintaining endogenous channel distribution (top two) or vice versa (bottom two) impacts optogenetic excitability. It is clear that the interaction of structural morphology and channel distribution has an impact on optogenetic excitability. The axon subregion appears to be the most susceptible with tripled excitation threshold of pyr₃ vs. pyr₂. Neuron degeneracy, i.e., the ability to perform the same functioning whilst being structurally different or having different ion channel distributions (Migliore et al., 2018), is thus something that should be taken into account in determining irradiance thresholds. This is also observed in the variability of I_th in Supplementary Figure 1. Still, the basal dendrites region is also the most excitable in pyr₃ and pyr₄. Combined with rank 1 and 2 for pyr₁ and pyr₂, respectively, it can be concluded that this is the most effective subcellular target region for opsin expression in CA1 pyramidal cells.

TABLE 5

Table 5. Mean relative difference in I_th due to three-dimensional structural morphology (top) or by changing channel distribution (bottom; relative with respect to original).

This spatial dependence was also observed in the study of Grossman et al. (2013). With the specific conductance (g_ChR2) as metric, they determined whole cell illumination to be most efficient, i.e., uniform opsin distribution over whole cell with a uniform light field compared to soma, axon initial segment or apical dendrite confinement. They determined that for a 20 ms pulse and irradiance of 1 mW/mm² the required g_ChR2 when in all sections is only 6% of that when restricted to the soma. On the other hand, G_max was 60% higher. These values are in agreement with our results where the ratio of the specific conductance of all sections to soma targeted expression is 2–4% under the same conditions in the pyramidal cells. However, in our study we argue that ranking should be based on G_max, i.e., where the number of opsin complexes is fixed. This translates toward an equal comparison of total elicited photocurrent, while, on the other hand, for a fixed g_ChR2 the total photocurrent is scaled by the surface area. As a result, the confinement to the soma is classified here to be more excitable.

After correction for the difference in rectification function [i.e., G(V = −68.83 mV) = 1 in Grossman et al. (2011) vs. 0.07 in this study], the absolute values of G_max were slightly lower but in the same order as reported in Grossman et al. (2013). For a G_max of 1 μS (= 0.07 μS after correction) and with a single channel conductance of 40–100 fS this translates toward expression of 0.71–1.77 10⁶ opsin complexes. Spread over the whole cell this is ±50 channels/μm² but confined only to the soma this rises toward >1,000 channels/μm². This value is higher than the estimated 130 channels/μm² based on bacteriorhodopsin expression (Nagel et al., 1995; Foutz et al., 2012) but lower than the indirectly estimated 4.4 10⁴ channels/μm² by Arlow et al. (2013). With our tested G_max values up to 100 μS especially when restricted to the soma, this could pose cellular toxicity problems, if these channel numbers would be achieved (Yizhar et al., 2011; Rost et al., 2017, 2022). This can be avoided if single channel conductance is increased.

The opsin is in all conditions uniformly distributed but can be restricted to a spatial region. In-vivo this highly specific and discrete separation is not possible. Still, by merging the opsin with signaling and targeting constructs, localized enhancement can be obtained. For instance, the addition of the soma-targeting motif of the soma- and proximal dendrite-localized voltage-gated potassium-channel Kv2.1 improves somato-dendritic expression (Rost et al., 2017; Mahn et al., 2018). The real distribution in those cases is not known. A normal distribution could be imposed but this would come with two more degrees of freedom. Therefore, the spatially restricted but uniform distribution is used. Our results encourage localized enhancement and advances in this research direction. A reduction of I_th with more than 64% can be achieved via subcellular specificity. This is tempered toward, but still significant, increases in SoFPAN of 33–100%, when light propagation is included. Consequently, if made possible, spatial confinement of opsins to specific membrane compartments could significantly increase optogenetic efficiency. On the other hand, more than 76% reduction of I_th between optimal and worst subcellular regions is possible. This is also reflected in SoFPAN, where an ideal subcellular-pitch combination can result in a 1.5–10 fold increase compared to the worst combination (cf. Table 3). A good knowledge of the optogenetic interaction at the subcellular level is therefore required in order to achieve the optimal configuration.

4.2. Optimal fiber position

Due to the finite size and discrete nature of the test grid (121 points), the SoFPAN, and optimal and worst positions could not be unambiguously determined. An upper and lower bound of the SoFPAN is defined (cf. Section 2.4.1), indicating a larger uncertainty at lower intensities. In case of the optimal fiber position, a tie-breaker based on the TAC is introduced. From a homeostatic point of view, it is ideal that the required result is achieved via the lowest perturbation of normal functioning. High/long transmembrane currents could lead to ion concentration imbalance. Especially when using ChR2(H134R), which has a high H⁺ permeability, this can result in neuron acidification, which in turn can result in decreased neuron functioning or unexpected behavior (Schneider et al., 2013; Mahn et al., 2016; Vierock et al., 2017). Therefore, the positions that generate the lowest TAC are preferred. The results in Figure 3C show that the optimal position is at the depth of the region of subcellular expression or with a focus on the most excitable region (in case of opsin distribution over the whole or majority of the cell). Overall, the rank of excitability between uniform and M.C is unchanged field stimulation. Subcellular excitability appears to be dominant over spatial distribution. This spatial preference was also observed in Foutz et al. (2012). In their L5 pyramidal model, they found the apical tuft and soma to be most excitable.

For the π/2 pitch, the optimal and worst position are the most stationary for increasing intensities. For the other pitches (0 and π), either the optimal or worst position is more smeared out while they are located more closely to each other for low intensities. Consequently, there is a higher risk for sub-optimal fiber positioning. Combined with the highest excitability according to the SoFPAN (cf. Figure 3E), it can be conclude that π/2 is the better fiber position.

4.3. Contribution of optical field simulation

This study combined simulations of light propagation and neuronal modeling. Light propagation is simulated using the Monte Carlo method for a uniform medium. The hippocampus is a predominantly gray matter structure. However, there is uncertainty on the exact values of the optical parameters (amplified by inter- and extra-polation). The effect of the uncertainty of these parameters is tested using the elementary effects method. The influence of ${\mu }_s^{\prime }$ on the excitability is ranked in the middle, while μ_a is ranked lower. On the optimal position they were ranked much higher. The median and maximum μ^* on SoFPAN are, respectively 0.14 and 1.71 mm² (0.04 and 1.19 mm²) for ${\mu }_s^{\prime }$ (μ_a). These parameters have some influence, but are subordinate to the other uncertain parameters such as subcellular location, expression level and cell morphology. Moreover, the need to include the light intensity profile in the neuron simulation was addressed by calculating the SoFPAN from excitation thresholds under uniform illumination, as well. Deviations of more than 20% are observed in more than 25% of the tested setups (soma excluded). Confinement to the basal dendrites result in the lowest percentage (18.17%) while the highest is achieved in case of the axon subcellular expression (34.60%). Overall, investigating optogenetic excitability under uniform field conditions provides a good initial approximation, but, accuracy drops for larger and asymmetrical section distributions. The latter was also observed in the strong pitch dependence of axon excitability in the pyramidal cells.

4.4. Limitations and future work

This study focuses on single pulse excitation of CA1 cells. The occurrence of other spiking behavior is excluded. Likewise, when calculating the SoFPAN with high I_fiber, a cell near the fiber end may exhibit bursting or depolarization block due to intense irradiance (Grossman et al., 2011, 2013). Unlike with electrical stimulation, the photocurrent saturates for high light intensities (see Figure 1C). Therefore, extreme behavior is not expected when short pulses are applied. Additionally, the studied ChR2(H134R) opsin exhibits light-dark adaptation, i.e., the photocurrent is higher for a full dark adapted neuron but decreases toward a steady state value under prolonged illumination. The recovery time is in the order of seconds. During pulsed stimulation, the photocurrent of the first pulse will be higher than that of the subsequent pulses. Therefore, higher irradiances will be required to reliably elicit pulse-locked spiking (Foutz et al., 2012). Future work should test the capability of eliciting reliable spiking when opsin expression is confined to a specific subcellular region. Additionally, since this study focuses on isolated cells that are at rest prior to optogenetic stimulation, it is necessary to investigate if phase locking is possible in a network setting, quantify its impact on excitation thresholds, observe subcellular excitation's effects on cellular and network responses (Youssef et al., 2023), and evaluate its influence on synaptic plasticity. Therefore, the interaction in neuronal networks will be of interest in future work, particularly with a focus on hyperexcitable systems such as those in temporal lobe epilepsy. Clearly, the cell's optogenetic excitability depends on multiple factors. The results in Table 5 show both 3D structural morphology and endogenous channel distribution dependencies. Still, the individual impact and interaction effects are yet to be determined. Obtaining a better understanding will be of interest in future work.

The effective level of opsin expression in-vivo is uncertain. Furthermore, to account for potential improvements in plasma membrane expression (Rost et al., 2017; Mahn et al., 2018), G_max is treated as a free parameter. However, due to the presence of inward rectification, it is unclear how the model parameter g_ChR2 relates to the actual opsin expression level in-vivo. In our formulation of G(V), all proportionality factors are absorbed by g_ChR2 (see Equations 1 and 2). While in the formulation of Grossman et al. (2011) and Williams et al. (2013), G(V = −68.83 mV) = 1 and G(V = −76.07 mV) = 1, respectively. These rectification functions cause a reduction of g_ChR2 with a factor of 14.15 or 12.89 at those specific membrane potentials, compared to our formulation in Equation 2.

As aforementioned, there is still some uncertainty on the tissue optical properties. Different studies have reported values that can differ up to an order of magnitude (Yaroslavsky et al., 2002; Gebhart et al., 2006; Johansson, 2010; Genina et al., 2019). Furthermore, brain tissue is binary classified as either gray or white matter, while tissue gradation is more continuous. The impact of the uncertainty of these parameters on the optogenetic excitability is tested in this study. However, it is investigated only locally near the parameters' reported means (see Table 1), while the reduced scattering coefficient of white matter is reported to be 7-times higher than that of gray matter (Yaroslavsky et al., 2002). Additionally, tissue alterations due to foreign body reaction occur. A fibrous capsule is formed around the implanted fiber as reaction to blood-brain-barrier injury and gliosis caused by the presence of the implanted fiber itself (Vandekerckhove et al., 2020). In future work light propagation in a heterogeneous medium and its effect on excitability could be determined. Mesh and voxel based Monte Carlo algorithms exist that can accurately compute the light field distribution in complex tissues (Yan and Fang, 2020). Moreover, modern deep learning algorithms can be utilized to reduce the inherent stochastic noise of these Monte Carlo simulations (Ardakani et al., 2022). Exploring the propagation of light at different wavelengths, such as for the excitation of red-shifted opsins, would also be of interest. The simulated fiber has a flat tip with fixed diameter and numerical aperture. Fiber tapering, flat tip patterning and alteration of the geometric properties can result in improved output coupling or broadened and multi-site illumination. However, this is out of the scope of this research (Kostovski et al., 2014; Vandekerckhove et al., 2020).

New opsins, either natural or genetically engineered, are discovered on a yearly basis (Yizhar et al., 2011; Klapoetke et al., 2014; Vierock et al., 2017, 2022; Oppermann et al., 2019). While this is generally beneficial, it can hinder the development of optogenetic tools. Dividing research among multiple opsins may limit the understanding of a single opsin's interactions with neurons and its capabilities. The question arises whether the gathered insights here are transferable to other opsins as well. Previous studies have demonstrated the influence of channel kinetics on factors like irradiance thresholds, spike reliability, and behavior (Grossman et al., 2011; Mattis et al., 2012; Klapoetke et al., 2014). The exact values of I_th will thus differ for another opsin. However, these values are already uncertain due to multiple other uncertainties in other parameters like G_max. Furthermore, these differences will affect all tested cases equally. Therefore, we expect that the observed trends and rankings regarding optogenetic excitability will be applicable across different opsins. In future work, a similar study could be performed focusing on optogenetic silencing with inhibitory opsins like GtACR2 and WiChr (Govorunova et al., 2015; Vierock et al., 2022).

Tissue illumination causes heating. To avoid permanent tissue damage, the local temperature increase cannot exceed 6°C (Acharya et al., 2021). Moreover, behavioral changes are already possible at lower temperature changes (1°C). Several neural parameters, e.g., capacitance, ion channel conductance, and transmitter release and uptake, have been shown to be temperature dependent (Kiyatkin, 2019; Fekete et al., 2020; Peixoto et al., 2020; Acharya et al., 2021). The reported SoFPAN values are for fiber intensities up to 1,000 mW/mm². After extrapolation of the change in temperature results reported in Acharya et al. (2021) (Figure 3A), this intensity corresponds with an estimated temperature increase of 4.85°C after 100 ms. Consequently, temperature-induced changes in optogenetic excitability should be included in future work.

5. Conclusion

In this in-silico study, we examined the optogenetic excitability of four CA1 cells using ChR2(H134R). Our findings reveal that, for a fixed amount of opsin channels (G_max), confining the opsin to specific neuronal membrane compartments significantly enhances excitability. This confinement leads to threshold reductions exceeding 64% and up to 100% gains in the surface of fiber positions for the activation of neurons. Additionally, we determined that the perpendicular orientation of the fiber relative to the somato-dendritic axis yields superior results. Furthermore, we observed substantial inter-cell variability, with differences in thresholds above 20%. The bistratified cell exhibited the least excitability, while pyramidal cell 1 demonstrated the highest excitability, especially when the opsin is confined to the basal dendrites. These findings highlight the importance of considering neuron degeneracy while developing optogenetic tools. By screening various uncertain parameters, we identified opsin location and G_max having the greatest impact on simulation outcomes. Our study showed the advantages of computational modeling coupled with light propagation. An increased excitability is seen with optimal fiber positioning, i.e., perpendicular to the somatic-dendritic axis and focus on the most excitable cell region. Spatial confinement and enhancements of opsin expression levels are promoted strategies to improve optogenetic excitability. However, it should be noted that uncertainty in these parameters limits determining the exact irradiance thresholds.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation. The used simulation tools are available on Github: https://github.com/Rschoete/3DOptogeneticStimulationCA1_public.

Author contributions

RS, TT, LW, WJ, RR, and ET contributed to the study design. RS performed all simulations. All authors aided in preparing the manuscript, read, and approved the final manuscript.

Funding

RS was a Ph.D. Fellow of the FWO-V (FR; Research Foundation Flanders, Belgium; Grant number: 11E4520N). TT was a postdoctoral fellow of the FWO-V (Grant number: 1230222N). This work was supported by the BOF project SOFTRESET.

Acknowledgments

The computational resources (Stevin Supercomputer Infrastructure) and services used in this work were provided by the VSC (Flemish Supercomputer Center), funded by Ghent University, FWO and the Flemish Government—department EWI.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fncom.2023.1229715/full#supplementary-material

References

Acharya, A. R., Vandekerckhove, B., Larsen, L. E., Delbeke, J., Wadman, W. J., Vonck, K., et al. (2021). In vivo blue light illumination for optogenetic inhibition: effect on local temperature and excitability of the rat hippocampus. J. Neural Eng. 18:066038. doi: 10.1088/1741-2552/ac3ef4

PubMed Abstract | CrossRef Full Text | Google Scholar

Ardakani, M. R., Yu, L., Kaeli, D. R., and Fang, Q. (2022). Framework for denoising monte carlo photon transport simulations using deep learning. J. Biomed. Optic. 27, 083019. doi: 10.1117/1.JBO.27.8.083019

PubMed Abstract | CrossRef Full Text | Google Scholar

Arlow, R., Foutz, T., and McIntyre, C. (2013). Theoretical principles underlying optical stimulation of myelinated axons expressing channelrhodopsin-2. Neuroscience 248, 541–551. doi: 10.1016/j.neuroscience.2013.06.031

PubMed Abstract | CrossRef Full Text | Google Scholar

Ascoli, G. A., Alonso-Nanclares, L., Anderson, S. A., Barrionuevo, G., Benavides-Piccione, R., Burkhalter, A., et al. (2008). Petilla terminology: nomenclature of features of gabaergic interneurons of the cerebral cortex. Nat. Rev. Neurosci. 9, 557–568. doi: 10.1038/nrn2402

PubMed Abstract | CrossRef Full Text | Google Scholar

Assaf, F., and Schiller, Y. (2016). The antiepileptic and ictogenic effects of optogenetic neurostimulation of PV-expressing interneurons. J. Neurophysiol. 116, 1694–1704. doi: 10.1152/jn.00744.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

Bansal, H., Gupta, N., and Roy, S. (2021). Theoretical analysis of optogenetic spiking with chrmine, breaches and CschRimson-expressing neurons for retinal prostheses. J. Neural Eng. 18:0460b8. doi: 10.1088/1741-2552/ac1175

PubMed Abstract | CrossRef Full Text | Google Scholar

Berndt, A., Schoenenberger, P., Mattis, J., Tye, K. M., Deisseroth, K., Hegemann, P., et al. (2011a). High-efficiency channelrhodopsins for fast neuronal stimulation at low light levels. Proc. Natl. Acad. Sci. U.S.A. 108, 7595–7600. doi: 10.1073/pnas.1017210108

PubMed Abstract | CrossRef Full Text | Google Scholar

Boyden, E. S., Zhang, F., Bamberg, E., Nagel, G., and Deisseroth, K. (2005). Millisecond-timescale, genetically targeted optical control of neural activity. Nat. Neurosci. 8, 1263–1268. doi: 10.1038/nn1525

PubMed Abstract | CrossRef Full Text | Google Scholar

Butson, C. R., Cooper, S. E., Henderson, J. M., and McIntyre, C. C. (2007). Patient-specific analysis of the volume of tissue activated during deep brain stimulation. NeuroImage 34, 661–670. doi: 10.1016/j.neuroimage.2006.09.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Campolongo, F., Saltelli, A., and Cariboni, J. (2011). From screening to quantitative sensitivity analysis. A unified approach. Comput. Phys. Commun. 182, 978–988. doi: 10.1016/j.cpc.2010.12.039

CrossRef Full Text | Google Scholar

Carnevale, N. T., and Hines, M. L. (2006). The NEURON Book. Cambridge: Cambridge University Press. doi: 10.1017/CBO9780511541612

CrossRef Full Text | Google Scholar

Carrette, S., Boon, P., Sprengers, M., Raedt, R., and Vonck, K. (2015). Responsive neurostimulation in epilepsy. Expert Rev. Neurotherap. 15, 1445–1454. doi: 10.1586/14737175.2015.1113875

PubMed Abstract | CrossRef Full Text | Google Scholar

Cela, E., McFarlan, A. R., Chung, A. J., Wang, T., Chierzi, S., Murai, K. K., et al. (2019). An optogenetic kindling model of neocortical epilepsy. Sci. Rep. 9:5236. doi: 10.1038/s41598-019-41533-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Cela, E., and Sjöström, P. J. (2019). Novel optogenetic approaches in epilepsy research. Front. Neurosci. 13:947. doi: 10.3389/fnins.2019.00947

PubMed Abstract | CrossRef Full Text | Google Scholar

Deisseroth, K. (2011). Optogenetics. Nat. Methods 8, 26–29. doi: 10.1038/nmeth.f.324

PubMed Abstract | CrossRef Full Text | Google Scholar

Duffley, G., Anderson, D. N., Vorwerk, J., Dorval, A. D., and Butson, C. R. (2019). Evaluation of methodologies for computing the deep brain stimulation volume of tissue activated. J. Neural Eng. 16:066024. doi: 10.1088/1741-2552/ab3c95

PubMed Abstract | CrossRef Full Text | Google Scholar

Fekete, Z., Horváth, A. C., and Zátonyi, A. (2020). Infrared neuromodulation:a neuroengineering perspective. J. Neural Eng. 17:051003. doi: 10.1088/1741-2552/abb3b2

PubMed Abstract | CrossRef Full Text | Google Scholar

Foutz, T. J., Arlow, R. L., and McIntyre, C. C. (2012). Theoretical principles underlying optical stimulation of a channelrhodopsin-2 positive pyramidal neuron. J. Neurophysiol. 107, 3235–3245. doi: 10.1152/jn.00501.2011

PubMed Abstract | CrossRef Full Text | Google Scholar

Gebhart, S. C., Lin, W. C., and Mahadevan-Jansen, A. (2006). In vitro determination of normal and neoplastic human brain tissue optical properties using inverse adding-doubling. Phys. Med. Biol. 51, 2011–2027. doi: 10.1088/0031-9155/51/8/004

PubMed Abstract | CrossRef Full Text | Google Scholar

Genina, E. A., Bashkatov, A. N., Tuchina, D. K., Timoshina, P. A. D., Navolokin, N., Shirokov, A., et al. (2019). Optical properties of brain tissues at the different stages of glioma development in rats: pilot study. Biomed. Opt. Exp. 10:5182. doi: 10.1364/BOE.10.005182

PubMed Abstract | CrossRef Full Text | Google Scholar

Govorunova, E. G., Sineshchekov, O. A., Janz, R., Liu, X., and Spudich, J. L. (2015). Natural light-gated anion channels: a family of microbial rhodopsins for advanced optogenetics. Science 349, 647–650. doi: 10.1126/science.aaa7484

PubMed Abstract | CrossRef Full Text | Google Scholar

Gradinaru, V., Thompson, K. R., and Deisseroth, K. (2008). ENPHR: a natronomonas halorhodopsin enhanced for optogenetic applications. Brain Cell Biol. 36, 129–139. doi: 10.1007/s11068-008-9027-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Gradinaru, V., Zhang, F., Ramakrishnan, C., Mattis, J., Prakash, R., Diester, I., et al. (2010). Molecular and cellular approaches for diversifying and extending optogenetics. Cell 141, 154–165. doi: 10.1016/j.cell.2010.02.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Grossman, N., Nikolic, K., Toumazou, C., and Degenaar, P. (2011). Modeling study of the light stimulation of a neuron cell with channelrhodopsin-2 mutants. IEEE Trans. Biomed. Eng. 58, 1742–1751. doi: 10.1109/TBME.2011.2114883

PubMed Abstract | CrossRef Full Text | Google Scholar

Grossman, N., Simiaki, V., Martinet, C., Toumazou, C., Schultz, S. R., and Nikolic, K. (2013). The spatial pattern of light determines the kinetics and modulates backpropagation of optogenetic action potentials. J. Comput. Neurosci. 34, 477–488. doi: 10.1007/s10827-012-0431-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Gupta, N., Bansal, H., and Roy, S. (2019). Theoretical optimization of high-frequency optogenetic spiking of red-shifted very fast-chrimson-expressing neurons. Neurophotonics 6:1. doi: 10.1117/1.NPh.6.2.025002

PubMed Abstract | CrossRef Full Text | Google Scholar

Johansson, J. D. (2010). Spectroscopic method for determination of the absorption coefficient in brain tissue. J. Biomed. Opt. 15:057005. doi: 10.1117/1.3495719

PubMed Abstract | CrossRef Full Text | Google Scholar

Kiyatkin, E. A. (2019). Brain temperature and its role in physiology and pathophysiology: lessons from 20 years of thermorecording. Temperature 6, 271–333. doi: 10.1080/23328940.2019.1691896

PubMed Abstract | CrossRef Full Text | Google Scholar

Klapoetke, N. C., Murata, Y., Kim, S. S., Pulver, S. R., Birdsey-Benson, A., Cho, Y. K., et al. (2014). Independent optical excitation of distinct neural populations. Nat. Methods 11, 338–346. doi: 10.1038/nmeth.2836

PubMed Abstract | CrossRef Full Text | Google Scholar

Kostovski, G., Stoddart, P. R., and Mitchell, A. (2014). The optical fiber tip: an inherently light-coupled microscopic platform for micro- and nanotechnologies. Adv. Mater. 26, 3798–3820. doi: 10.1002/adma.201304605

PubMed Abstract | CrossRef Full Text | Google Scholar

Krook-Magnuson, E., Armstrong, C., Bui, A., Lew, S., Oijala, M., and Soltesz, I. (2015). In vivo evaluation of the dentate gate theory in epilepsy. J. Physiol. 593, 2379–2388. doi: 10.1113/JP270056

PubMed Abstract | CrossRef Full Text | Google Scholar

Krook-Magnuson, E., Armstrong, C., Oijala, M., and Soltesz, I. (2013). On-demand optogenetic control of spontaneous seizures in temporal lobe epilepsy. Nat. Commun. 4:1376. doi: 10.1038/ncomms2376

PubMed Abstract | CrossRef Full Text | Google Scholar

Ladas, T. P., Chiang, C.-C., Gonzalez-Reyes, L. E., Nowak, T., and Durand, D. M. (2015). Seizure reduction through interneuron-mediated entrainment using low frequency optical stimulation. Exp. Neurol. 269, 120–132. doi: 10.1016/j.expneurol.2015.04.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Lehtinen, K., Nokia, M. S., and Takala, H. (2022). Red light optogenetics in neuroscience. Front. Cell. Neurosci. 15:532. doi: 10.3389/fncel.2021.778900

PubMed Abstract | CrossRef Full Text | Google Scholar

Lu, Y., Zhong, C., Wang, L., Wei, P., He, W., Huang, K., et al. (2016). Optogenetic dissection of ictal propagation in the hippocampal-entorhinal cortex structures. Nat. Commun. 7:10962. doi: 10.1038/ncomms10962

PubMed Abstract | CrossRef Full Text | Google Scholar

Mahn, M., Gibor, L., Patil, P., Malina, K. C.-K., Oring, S., Printz, Y., et al. (2018). High-efficiency optogenetic silencing with soma-targeted anion-conducting channelrhodopsins. Nat. Commun. 9:4125. doi: 10.1038/s41467-018-06511-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Mahn, M., Prigge, M., Ron, S., Levy, R., and Yizhar, O. (2016). Biophysical constraints of optogenetic inhibition at presynaptic terminals. Nat. Neurosci. 19, 554–556. doi: 10.1038/nn.4266

PubMed Abstract | CrossRef Full Text | Google Scholar

Markram, H., Muller, E., Ramaswamy, S., Reimann, M., Abdellah, M., Sanchez, C., et al. (2015). Reconstruction and simulation of neocortical microcircuitry. Cell 163, 456–492. doi: 10.1016/j.cell.2015.09.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Mattis, J., Tye, K. M., Ferenczi, E. A., Ramakrishnan, C., O'Shea, D. J., Prakash, R., et al. (2012). Principles for applying optogenetic tools derived from direct comparative analysis of microbial opsins. Nat. Methods 9, 159–172. doi: 10.1038/nmeth.1808

PubMed Abstract | CrossRef Full Text | Google Scholar

Messier, J. E., Chen, H., Cai, Z.-L., and Xue, M. (2018). Targeting light-gated chloride channels to neuronal somatodendritic domain reduces their excitatory effect in the axon. eLife 7:e38506. doi: 10.7554/eLife.38506.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Migliore, R., Lupascu, C. A., Bologna, L. L., Romani, A., Courcol, J.-D., Antonel, S., et al. (2018). The physiological variability of channel density in hippocampal ca1 pyramidal cells and interneurons explored using a unified data-driven modeling workflow. PLoS Comput. Biol. 14:e1006423. doi: 10.1371/journal.pcbi.1006423

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagel, G., Möckel, B., Büldt, G., and Bamberg, E. (1995). Functional expression of bacteriorhodopsin in oocytes allows direct measurement of voltage dependence of light induced h+ pumping. FEBS Lett. 377, 263–266. doi: 10.1016/0014-5793(95)01356-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Nikolic, K., Degenaar, P., and Toumazou, C. (2006).“Modeling and engineering aspects of channelrhodopsin2 system for neural photostimulation,” in 2006 International Conference of the IEEE Engineering in Medicine and Biology Society (New York, NY: IEEE), 1626–1629. doi: 10.1109/IEMBS.2006.260766

PubMed Abstract | CrossRef Full Text | Google Scholar

Nikolic, K., Grossman, N., Grubb, M. S., Burrone, J., Toumazou, C., and Degenaar, P. (2009). Photocycles of channelrhodopsin-2. Photochem. Photobiol. 85, 400–411. doi: 10.1111/j.1751-1097.2008.00460.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Oppermann, J., Fischer, P., Silapetere, A., Liepe, B., Rodriguez-Rozada, S., Flores-Uribe, J., et al. (2019). Mermaids: a family of metagenomically discovered marine anion-conducting and intensely desensitizing channelrhodopsins. Nat. Commun. 10:3315. doi: 10.1038/s41467-019-11322-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Owen, S. F., Liu, M. H., and Kreitzer, A. C. (2019). Thermal constraints on in vivo optogenetic manipulations. Nat. Neurosci. 22, 1061–1065. doi: 10.1038/s41593-019-0422-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Peixoto, H. M., Cruz, R. M. S., Moulin, T. C., and Leão, R. N. (2020). Modeling the effect of temperature on membrane response of light stimulation in optogenetically-targeted neurons. Front. Comput. Neurosci. 14:5. doi: 10.3389/fncom.2020.00005

PubMed Abstract | CrossRef Full Text | Google Scholar

Reilly, J. P. (1998). “Excitation models,” in Applied Bioelectricity (New York, NY: Springer). doi: 10.1007/978-1-4612-1664-3_4

CrossRef Full Text | Google Scholar

Romani, A., Schürmann, F., Markram, H., and Migliore, M. (2022). “Reconstruction of the hippocampus,” in Computational Modelling of the Brain. Advances in Experimental Medicine and Biology, Vol. 1359, eds M. Giugliano, M. Negrello, and D. Linaro (Cham: Springer). doi: 10.1007/978-3-030-89439-9_11

PubMed Abstract | CrossRef Full Text | Google Scholar

Rost, B. R., Schneider-Warme, F., Schmitz, D., and Hegemann, P. (2017). Optogenetic tools for subcellular applications in neuroscience. Neuron 96, 572–603. doi: 10.1016/j.neuron.2017.09.047

PubMed Abstract | CrossRef Full Text | Google Scholar

Rost, B. R., Wietek, J., Yizhar, O., and Schmitz, D. (2022). Optogenetics at the presynapse. Nat. Neurosci. 25, 984–998. doi: 10.1038/s41593-022-01113-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Saltelli, A., Aleksankina, K., Becker, W., Fennell, P., Ferretti, F., Holst, N., et al. (2019). Why so many published sensitivity analyses are false: a systematic review of sensitivity analysis practices. Environ. Model. Softw. 114, 29–39. doi: 10.1016/j.envsoft.2019.01.012

CrossRef Full Text | Google Scholar

Saray, S., Rossert, C. A., Appukuttan, S., Migliore, R., Vitale, P., Lupascu, C. A., et al. (2021). HippoUnit: a software tool for the automated testing and systematic comparison of detailed models of hippocampal neurons based on electrophysiological data. PLoS Comput. Biol. 17:e1008114. doi: 10.1101/2020.07.02.184333

PubMed Abstract | CrossRef Full Text | Google Scholar

Schneider, F., Gradmann, D., and Hegemann, P. (2013). Ion selectivity and competition in channelrhodopsins. Biophys. J. 105, 91–100. doi: 10.1016/j.bpj.2013.05.042

PubMed Abstract | CrossRef Full Text | Google Scholar

Schoeters, R., Tarnaud, T., Martens, L., Joseph, W., Raedt, R., and Tanghe, E. (2021). Double two-state opsin model with autonomous parameter inference. Front. Comput. Neurosci. 15:688331. doi: 10.3389/fncom.2021.688331

PubMed Abstract | CrossRef Full Text | Google Scholar

Shin, Y., and Kwon, H.-S. (2016). Mesh-based Monte Carlo method for fibre-optic optogenetic neural stimulation with direct photon flux recording strategy. Phys. Med. Biol. 61, 2265–2282. doi: 10.1088/0031-9155/61/6/2265

PubMed Abstract | CrossRef Full Text | Google Scholar

Sprengers, M., Raedt, R., Larsen, L. E., Delbeke, J., Wadman, W. J., Boon, P., et al. (2020). Deep brain stimulation reduces evoked potentials with a dual time course in freely moving rats: potential neurophysiological basis for intermittent as an alternative to continuous stimulation. Epilepsia 61, 903–913. doi: 10.1111/epi.16498

PubMed Abstract | CrossRef Full Text | Google Scholar

Stujenske, J., Spellman, T., and Gordon, J. (2015). Modeling the spatiotemporal dynamics of light and heat propagation for in-vivo optogenetics. Cell Rep. 12, 525–534. doi: 10.1016/j.celrep.2015.06.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Tønnesen, J., and Kokaia, M. (2017). Epilepsy and optogenetics: can seizures be controlled by light? Clin. Sci. 131, 1605–1616. doi: 10.1042/CS20160492

PubMed Abstract | CrossRef Full Text | Google Scholar

Tønnesen, J., Sørensen, A. T., Deisseroth, K., Lundberg, C., and Kokaia, M. (2009). Optogenetic control of epileptiform activity. Proc. Natl. Acad. Sci. U.S.A. 106, 12162–12167. doi: 10.1073/pnas.0901915106

PubMed Abstract | CrossRef Full Text | Google Scholar

Vandekerckhove, B., Missinne, J., Vonck, K., Bauwens, P., Verplancke, R., Boon, P., et al. (2020). Technological challenges in the development of optogenetic closed-loop therapy approaches in epilepsy and related network disorders of the brain. Micromachines 12:38. doi: 10.3390/mi12010038

PubMed Abstract | CrossRef Full Text | Google Scholar

Vierock, J., Grimm, C., Nitzan, N., and Hegemann, P. (2017). Molecular determinants of proton selectivity and gating in the red-light activated channelrhodopsin chrimson. Sci. Rep. 7:9928. doi: 10.1038/s41598-017-09600-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Vierock, J., Peter, E., Grimm, C., Rozenberg, A., Chen, I.-W., Tillert, L., et al. (2022). Wichr, a highly potassium-selective channelrhodopsin for low-light one- and two-photon inhibition of excitable cells. Sci. Adv. 8:eadd7729. doi: 10.1126/sciadv.add7729

PubMed Abstract | CrossRef Full Text | Google Scholar

Walker, M. C., and Kullmann, D. M. (2020). Optogenetic and chemogenetic therapies for epilepsy. Neuropharmacology 168:107751. doi: 10.1016/j.neuropharm.2019.107751

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, L., Jacques, S. L., and Zheng, L. (1995). MCML-Monte Carlo modeling of light transport in multi-layered tissues. Comput. Methods Prog. Biomed. 47, 131–146. doi: 10.1016/0169-2607(95)01640-F

PubMed Abstract | CrossRef Full Text | Google Scholar

White, M., Mackay, M., and Whittaker, R. G. (2020). Taking optogenetics into the human brain: Opportunities and challenges in clinical trial design. Open Access J. Clin. Trials 12, 33–41. doi: 10.2147/OAJCT.S259702

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, J., and Entcheva, E. (2015). Optogenetic versus electrical stimulation of human cardiomyocytes: modeling insights. Biophys. J. 108, 1934–1945. doi: 10.1016/j.bpj.2015.03.032

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, J. C., Xu, J., Lu, Z., Klimas, A., Chen, X., Ambrosi, C. M., et al. (2013). Computational optogenetics: empirically-derived voltage- and light-sensitive channelrhodopsin-2 model. PLoS Comput. Biol. 9:e1003220. doi: 10.1371/journal.pcbi.1003220

PubMed Abstract | CrossRef Full Text | Google Scholar

Wykes, R. C., Kullmann, D. M., Pavlov, I., and Magloire, V. (2016). Optogenetic approaches to treat epilepsy. J. Neurosci. Methods 260, 215–220. doi: 10.1016/j.jneumeth.2015.06.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, S., and Fang, Q. (2020). Hybrid mesh and voxel based Monte Carlo algorithm for accurate and efficient photon transport modeling in complex bio-tissues. Biomed. Opt. Express 11:6262. doi: 10.1364/BOE.409468

PubMed Abstract | CrossRef Full Text | Google Scholar

Yaroslavsky, A. N., Schulze, P. C., Yaroslavsky, I. V., Schober, R., Ulrich, F., and Schwarzmaier, H.-J. (2002). Optical properties of selected native and coagulated human brain tissues in vitro in the visible and near infrared spectral range. Phys. Med. Biol. 47:305. doi: 10.1088/0031-9155/47/12/305

PubMed Abstract | CrossRef Full Text | Google Scholar

Yizhar, O., Fenno, L., Davidson, T., Mogri, M., and Deisseroth, K. (2011). Optogenetics in neural systems. Neuron 71, 9–34. doi: 10.1016/j.neuron.2011.06.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Youssef, D., Wittig, J. H., Jackson, S., Inati, S. K., and Zaghloul, K. A. (2023). Neuronal spiking responses to direct electrical microstimulation in the human cortex. J. Neurosci. 43, 4448–4460. doi: 10.1523/JNEUROSCI.1666-22.2023

PubMed Abstract | CrossRef Full Text | Google Scholar