**Changes in Ion Selectivity Following the Asymmetrical Addition of Charge to the Selectivity Filter of Bacterial Sodium Channels**

#### **Olena A. Fedorenko 1,\* ,**† **, Igor A. Khovanov <sup>2</sup> , Stephen K. Roberts <sup>1</sup> and Carlo Guardiani 3,**‡


Received: 28 October 2020; Accepted: 7 December 2020; Published: 9 December 2020

**Abstract:** Voltage-gated sodium channels (NaVs) play fundamental roles in eukaryotes, but their exceptional size hinders their structural resolution. Bacterial NaVs are simplified homologues of their eukaryotic counterparts, but their use as models of eukaryotic Na<sup>+</sup> channels is limited by their homotetrameric structure at odds with the asymmetric Selectivity Filter (SF) of eukaryotic NaVs. This work aims at mimicking the SF of eukaryotic NaVs by engineering radial asymmetry into the SF of bacterial channels. This goal was pursued with two approaches: the co-expression of different monomers of the NaChBac bacterial channel to induce the random assembly of heterotetramers, and the concatenation of four bacterial monomers to form a concatemer that can be targeted by site-specific mutagenesis. Patch-clamp measurements and Molecular Dynamics simulations showed that an additional gating charge in the SF leads to a significant increase in Na<sup>+</sup> and a modest increase in the Ca2<sup>+</sup> conductance in the NavMs concatemer in agreement with the behavior of the population of random heterotetramers with the highest proportion of channels with charge −5*e*. We thus showed that charge, despite being important, is not the only determinant of conduction and selectivity, and we created new tools extending the use of bacterial channels as models of eukaryotic counterparts.

**Keywords:** ion channel; selectivity; permeability; patch-clamp; computer simulations

#### **1. Introduction**

Voltage-gated sodium and calcium channels (NaVs and CaVs, respectively) are involved in a multitude of processes, including electrical signaling, secretion, and synaptic transmission [1]. The malfunction or dysregulation of NaVs and CaVs leads to a wide range of neurological, cardiovascular, and muscular disorders, including periodic paralysis [2], arrhythmia [3], and epilepsy [4], which highlights the importance of these molecules.

Eukaryotic NaVs and CaVs have similar structure and comprise a pore-forming α1 subunit of approximately 190–250 kDa, which co-assembles with a number of auxiliary subunits. The α1 subunit is organized in four domains, each comprising a voltage sensor (encompassing helices S1–S4) and a pore domain (including helices S5–S6) [5–9]. The four domains arranged around the pore are not identical, resulting in a channel structure that is asymmetric and pseudo-tetrameric [10]. The atomic level resolution of the structure of these molecules is essential to understand their structure-function relationships, but this task is particularly challenging for eukaryotic NaV channels, due to them being membrane integral proteins [11] and is exacerbated by their particularly large size. Consequently, to date only the structure of a single eukaryotic NaV has been resolved at the atomic level (3.8 Å), but the channel was not electrophysiologically characterized [12].

Bacterial NaVs and CaVs are simplified homologues of eukaryotic channels. Theyarehomotetrameric channels formed by four identical monomers corresponding to the four domains of the α1 subunit of eukaryotic channels [13–16]. Their minimalist structure has enabled the determination of high-resolution atomistic structures, which has allowed extensive structure-functional characterization with respect to their cation selectivity, gating, and binding of anesthetics e.g., [17–19]. Although a complete understanding of selectivity has not been arrived at, the availability of high-resolution structures has provided a detailed understanding of the atomistic-level interaction of ions in the Selectivity Filter (SF). Despite this, there is still no widely accepted predictive model for cation selectivity, representing gaps in our knowledge of the molecular mechanism of ion permeation and selectivity. Mutation studies suggest the fixed charge (Q*<sup>f</sup>* ) of the SF to be one of the major determinants of selectivity and permeation [15,20–32]. The Q*<sup>f</sup>* charge is at the core of many theoretical models attempting to explain the physical origins of cation selectivity in ion channels. For example, cation conduction has recently been modelled within the framework of Ionic Coulomb Blockade (ICB) [24,25], an electrostatic model with the aim of predicting Na<sup>+</sup> and Ca2<sup>+</sup> permeability based on knowing the actual Q*<sup>f</sup>* value of the SF.

Although prokaryotic and eukaryotic channels show the same general architecture along the axis of the pore (an outer vestibule and an inner water filled cavity separated by a narrow SF), the use of prokaryotic channels as models of their eukaryotic counterparts is limited by their lack of radial asymmetry. In the case of the homotetrameric bacterial NaVs, four identical monomers form the channel pore; in contrast, their eukaryotic counterparts are composed of four non-identical domains which introduce significant radial asymmetry [14,26–31]. This difference becomes evident at the level of the SF, where conduction and selectivity are controlled by a DEKA ring (Q*<sup>f</sup>* = −1*e*) in eukaryotic NaVs and an EEEE ring (Q*<sup>f</sup>* = −4*e*) in prokaryotic NaVs. Another puzzling fact is that the EEEE locus is typical of bacterial sodium-selective channels but also characterizes calcium-selective eukaryotic channels [16,32–34], thus leading to NaChBac being initially predicted to be Ca2+-selective. The existence of disparate sequences indicates that bacterial and eukaryotic channels enforce their ion preferences through different molecular strategies [15,34,35]. As a result, the selectivity and conduction mechanisms discovered in prokaryotes are not readily transferable to eukaryotes.

The puzzling functional similarity between bacterial NaVs and eukaryotic CaVs has been termed the "EEEE paradox" [36]. The paradox arises as a result of the violation of the assumption that Q*<sup>f</sup>* is the main driving force of cation selectivity. A possible resolution of the paradox is related to the existence of a conserved D residue in domain 2 of CaVs in the neighborhood of the EEEE ring. Monte Carlo simulations predicted this D residue (termed D2p51 in [37]) to occupy a position in close proximity to the EEEE locus. This observation led to the hypothesis that the locus imparting Ca2<sup>+</sup> permeability is actually EEEED, with a Q*<sup>f</sup>* value of −5*e* [25,37]. Moreover, when this conserved D residue in domain 2 of Cav1.2 (referred to as D707 in [29]) was replaced with neutral residues, a striking reduction in Ca2<sup>+</sup> binding to the SF was measured [29]. These results suggest D707 to be an important cation binding determinant of eukaryotic channels.

The different behavior of prokaryotic and eukaryotic voltage-gated sodium and calcium channels highlights the importance of incorporating radial asymmetry in the SF of prokaryotic channels. In our previous work [38] we reported the creation of a concatenated bacterial NaV, in which four NavMs monomers were covalently linked to form a stable single polypeptide chain, resembling the general structure of a eukaryotic NaV. This allowed the targeted mutagenesis of individual domains introducing radial asymmetry in the bacterial channel with the aim to gain further insight on the role of Q*<sup>f</sup>* as a determinant of ion selectivity. In the present study, we report the first attempt to mutate the concatemer and generate a bacterial sodium channel with radial asymmetry in the SF. In order to obtain atomistic-level detail of selectivity and permeation, the electrophysiological characterization

was integrated with Molecular Dynamics (MD) simulations of wild-type NavMs (Q*<sup>f</sup>* = −4*e*) and a mutant with an additional negative charge in the SF (Q*<sup>f</sup>* = −5*e*).

In the present study, we have also employed an independent yet complementary approach to introduce radial asymmetry into the SF of a bacterial sodium channel. Namely, a number of combinations of NaChBac monomers (differing in their amino acid composition and Q*<sup>f</sup>* value of the SF) were transfected into Chinese Hamster Ovary (CHO) cells to generate a random population of heterotetrameric channels with radial asymmetry in the SF. The mixed monomer approach using NaChBac monomers showed that Ca2<sup>+</sup> conduction is increased in channels with a Q*<sup>f</sup>* <sup>&</sup>gt; <sup>−</sup>4*<sup>e</sup>* (consistent with the proposed explanation for the EEEE paradox). Our data confirm the key role of the SF charge as the major determinant of conduction and selectivity. However, the failure to completely overturn the sodium selectivity of the NavMs concatemer to Ca2<sup>+</sup> selectivity (with much smaller relative Ca2<sup>+</sup> permeability exhibited by the −5*e* mutant NavMs concatemers compared to that for eukaryotic CaVs) suggests the existence of fine-tuning mechanisms of structural origin.

#### **2. Materials and Methods**

#### *2.1. Materials Generation of Mutant Bacterial Channels*

cDNA constructs encoding NaChBac (GenBank accession number BAB05220) and NavMs (GenBank accession number WP\_011712479) bacterial sodium channels were synthesized by EPOCH Life Science (www.epochlifescience.com). NavMs concatemer was subcloned into pTRACER-CMV2 (Invitrogen) downstream of CMV promoter, as described previously [38].

Site-directed mutagenesis was performed using specific primers containing the sequence for the desired amino acid substitutions (according to Q5® Site-Directed Mutagenesis Kit; New England BioLabs Inc., Hitchin, UK). For the generation of LEDWAS mutant from wild-type NaChBac, we used the forward primer CACGCTAGAGgatTGGGCGAGCG and the reversed primer ACCACTTGGAACAATGTTAAC; for LASWAS mutant, we used the forward primer GGTCACGCTAgccTCATGGGCGAGcggc and the reversed primer ACTTGGAACAATGTTAACAAACtaagc.

Q*<sup>f</sup>* = −5*e* NavMs mutants were generated from NavMs concatemer, which was designed with restriction sites delimiting each domain (Supplementary Figure S3A). Domain I (KpnI/EcoRI) and Domain II (EcoRI/EcoRV) were excised by restriction digest. The domain fragments were re-amplified by PCR using primer pairs to regenerate the restriction site prior to subcloning into vector pCR Blunt II-TOPO (Invitrogen): primers are Kpn1\_NavMs\_F (CCCGGTACCAGCCGCCA CCATGTCACGCAAAATAAG)/EcoRI\_ NavMs\_R (CCCGAATTCGGGCTCGTCCTCCCAGATG) for Domain I and EcoRI\_ NavMs\_F (CCCGAATTCATGTCTAGGAAGATCC)/EcoRV\_ NavMs\_F (CCCGATATCGGGCTCGTCCTCCCAGATG) for Domain II. Site-directed mutagenesis (for S179D, according to NavMs monomer residue nomenclature) was performed on each domain using primers LEDWSM\_ NavMs\_F (GACCTTAGAGgatTGGTCTATGGGC) and LEDWSM\_ NavMs\_R (ATCACCTGAAATAGTGTG) prior to the restriction enzyme-mediated excision and ligation (T4 DNA ligase; NEB) of the Domain DNA fragment in the NavMs concatemer at sites KpnI/EcoRI (for Domain I) and EcoRI/EcoRV (Domain II).

All the clones were sequenced to check for correct construction and to ensure that no unwanted PCR induced mutations had been introduced. DNA for the transfection of cells was prepared using Midi Plasmid Kit (Qiagen, Manchester, UK).

#### *2.2. Cell Culture and Transfection*

Chinese hamster ovary (CHO) and human embryonic kidney (HEK293T) cell lines were obtained from Dr. Stephen K. Roberts. Cells were cultured in DMEM high glucose with L-glutamine (Lonza, Slough, UK ) supplemented with 10% Fetal Bovine Serum (Thermo Scientific, Loughborough, UK) with the addition of 50 U/mL of penicillin and 50 µg/mL of streptomycin (Sigma, Irvine, UK). Cells were

maintained in a T25 flask (Thermo Scientific) at 37 ◦C in a 5% CO<sup>2</sup> incubator and passaged twice a week. Then, 24 h before transfection, the cells were seeded in 6-well plates (Corning, Deeside, UK) containing No.1 coverslips (Scientific Laboratory Supplies, Nottingham, UK). A total of 10 µL of transfection reagent (Mirus, Cardiff, UK) and 5 µg of plasmid DNA or a mixture of DNAs in defined proportions were equilibrated separately in 250 µL of UltraMEM™ Reduced Serum Medium (Lonza, Slough, UK) at room temperature for 5 min and then mixed and incubated at room temperature for 20 min to form the DNA-reagent complex. Treated cells (at 80% confluency) were supplemented with DNA-reagent complex and incubated at 37 ◦C and 5% CO<sup>2</sup> for 24–48 h before experiments.

#### *2.3. Electrophysiology*

Whole-cell voltage clamp recordings were performed at room temperature (20 ◦C) using an Axopatch 200A (Molecular Devices, Inc., Wokingham, UK) amplifier. Patch-clamp pipettes were pulled from borosilicate glass (Kimax, Kimble Company, Dover, USA) to resistances of 2–3 MOhm. Shanks of the pipette's tip were coated with bee's wax to reduce the pipette capacitance. The pipette solution contained (in mM) 15 Na-gluconate, 5 NaCl, 90 NMDG, 10 EGTA, and 20 HEPES, pH 7.4 adjustedwith 3 mM HCl). To record Na<sup>+</sup> influx currents the bath solution was (in mM) 140 Na-methanesulfonate, 5 CsCl, 10 HEPES and 10 glucose (pH 7.4 adjusted with 4.8 mM CsOH); for the measurement of the Ca2<sup>+</sup> influx currents, 140 mM of Na-methanesulfonate was replaced with 100 mM of Ca-methanesulfonate.

Data collection was initiated 3 min after obtaining the whole cell configuration to ensure the complete equilibration of the pipette solution and cytosol. The bath solution was grounded using a 3 M KCl agar bridge; the liquid junction potential determined experimentally (as described by [39]) agreed with that calculated (using JPCalc program, Clampex, Axon Instruments, Inc., Wokingham, UK) and was negligible. To ensure the complete exchange of the bath solution, electrophysiological recordings were initiated after >4 min of solution change. The rate of the gravity-fed perfusion system for the bath solution exchange was approximately 0.7 mL/min in a chamber volume of approximately 200 µL.

The results were analyzed using the Clampfit 10.1 software (Molecular Devices, Wokingham, UK) and OriginPro8 (OriginLab Corporation, Wellesley, MA, USA). Pooled data are presented as means ± SEM (*n*), where *n* is the number of independent experiments.

#### *2.4. Equilibrium Simulations of NavMs Channel*

The initial structure of wild-type NavMs was taken from the Protein Data Bank (PDB ID: 3ZJZ). Mutation S179D on chain A, and embedding in a membrane of 248 POPC molecules was performed using the CHARMM membrane builder [40,41]. The membrane was bathed on both sides by a 0.14 M NaCl solution or a 0.1 M CaCl<sup>2</sup> solution. The size of the simulation box was 102 × 102 × 86 Å and the total number of atoms in the four simulated systems was a little short of 90,000. All the acidic residues have been assigned a charge −1*e*, while basic residues have been assigned a charge +1*e* based on an analysis of the pKa values with the PROPKA program (server.poissonboltzmann.org/pdb2pqr). All the simulations were performed with the NAMD 2.11b2 [42] suite of programs using the ff14SB [43] force field for the protein and the Lipid14 force-field [44] for the phospholipids. As already observed in [45], in the absence of harmonic restraints the pore rapidly closes at the cytoplasmic gate. In order to avoid this behavior that likely results from the absence of the Voltage Sensor Domain in the simulated system, harmonic restraints (50 kcal/mol/Å<sup>2</sup> ) were applied to the backbone atoms of the transmembrane helices (residues 131–154 and 194–222) throughout the simulation. The four systems first underwent 10,000 steps of conjugate gradient minimization.

During equilibration, harmonic restraints were applied to non-hydrogen atoms of the protein backbone and side-chains (outside the transmembrane helices; residues 155–193), as well as to the phospholipid heads. A harmonic restraint was also applied to the dihedral angle formed by carbons 8, 9, 10, 11 of oleoyl acid and to the improper dihedral C1–C3–C2–O2 involving the three carbons of the glycerol unit and the hydroxyl oxygen linked to its central carbon. The equilibration was organized in six stages, whereby the constraints were gradually released. The values of the force constants used in the six stages can be found in Supplementary Table S1. The production run was carried out in the isothermal isobaric (NPT) ensemble for 100 ns (in NaCl solution) or 150 ns (in CaCl<sup>2</sup> solution). The pressure was kept at 1 atm by the Nose–Hoover Langevin piston method, while the temperature was kept at 300 K by coupling to a Langevin thermostat with a damping coefficient of 1 ps−<sup>1</sup> . Long-range electrostatic interactions were evaluated with the smooth particle mesh Ewald algorithm. For the short-range non-bonded interactions, we used a cutoff of 12 Å with a switching function at 10.0 Å. The integration time step was 2 fs, and the bonds between hydrogen and heavy atoms were fixed to eliminate the most rapid oscillatory motions. The Potential of Mean Force (PMF) was computed using equation *F*(*z*) = −*kBT* log(ρ(*z*)/ρ*<sup>b</sup>* ), where *k<sup>B</sup>* is the Boltzmann constant, *T* is the absolute temperature, ρ(*z*) is the density profile of sodium or calcium ions, and ρ*<sup>b</sup>* is the density of these ions in the bulk. Since the ion density in the channel is typically higher than in the bulk, the PMF normally has negative values. To avoid a divergence in the logarithmic expression of the PMF, we assigned *F*(*z*) = 0 when ρ(*z*) = 0—that is, in the regions of the channel that are never visited by ions.

#### *2.5. Current-Voltage Curves Calculation*

Current–voltage (IV) curves in NavMs were attained using the collective diffusion model introduced in [46], where the time-course *Q*(*t*) of the net charge transported across the channel at equilibrium is thought of as an unbiased random walk. The net charge transported in the time interval ∆*t* between two consecutive frames of the trajectory is ∆*Q* = P *z*1≤*z*≤*z*<sup>2</sup> *ei*∆*z<sup>i</sup> Lz* , where the sum runs over all ions *i*, such that *z*<sup>1</sup> ≤ *z* ≤ *z*2, *z*<sup>1</sup> = −4.5 Å and *z*<sup>2</sup> = 16.5 Å are the axial limits of the filter region somewhat extended in the vestibule and central cavity. The use of this extended SF gives us the opportunity to exploit the fluctuations due to ions exploring the vestibule region without entering into the SF as well as the aborted permeation events where the ion crosses the mouth of the SF but is immediately pulled back in due to the attraction of the acidic residues. In the expression, ∆*z<sup>i</sup>* is the axial displacement of the ion in the time interval ∆*t* and *L<sup>z</sup>* = *z*<sup>2</sup> − *z*<sup>1</sup> is the length of the SF. The time course of the charge, *Q*(*t*), can then be attained as *Q*(*t*) = P *ti*<*t* ∆*Q*(*ti*).

Diffusion theory predicts that, for sufficiently long times, the mean square displacement of the charge D *Q*2 (*t*) E grows linearly with a slope proportional to the diffusion coefficient, D *Q*2 (*t*) E ∼ 2*DQt* + *Const*. Applying linear response theory, the steady current induced by a small constant voltage *V* can be computed as *Isteady* = *DQV*/*kBT*. Using such an approach, the linear region of an IV curve can be computed based on the spontaneous ion fluctuations at equilibrium in the absence of any applied electric field.

#### **3. Results**

#### *3.1. Experimental Results*

To introduce radial asymmetry in the SF of NaChBac, two approaches were adopted. First, mixed populations of NaChBac monomers (differing in their amino acid composition and the Q*<sup>f</sup>* value of the SF) were co-transfected into CHO cells to generate hetero-tetrameric channels exhibiting radial asymmetry in the SFs. Second, we used a concatenated NavMs tetramer [38] to generate radial asymmetry in the SF by the targeted mutation of one of the four repeats.

### 3.1.1. Na+/Ca2<sup>+</sup> Selectivity for Randomly Mixed Populations of NaChBac Monomers

The random assembly of channel tetramers can be demonstrated taking advantage of the different electrophysiological properties of WT NaChBac and the L226P mutant illustrated by the recordings in Supplementary Figure S1. The L226P mutation causes conspicuous alterations in channel gating of NaChBac from depolarization-activated whole-cell currents to non-inactivating hyperpolarization-activated whole-cell currents [47], (Supplementary Figure S1A,B). The mutation shifts the voltage of the maximal current from −10 mV in WT to −180 mV in the mutant (Supplementary Figure S1D), thus currents at −10 mV originating from separate channel populations of WT and L226P homotetramers can be easily separated. The current recordings from CHO cells co-transfected with NaChBac-encoding WT:L226P cDNAs in a ratio 3:1 (Supplementary Figure S1C,E) exhibited unique currents at −10 mV, which cannot be explained by the simple addition of whole current traces from homotetramer channels formed from either L226P or wild-type NaChBac (note that there is no current at −10 mV from L226P channels), indicating that unique heterotetramers are being formed. Assuming that the assembly of heterotetramers is formed without bias, the proportions of channel types can be determined by binomial distribution. It is noteworthy that this assumption is in agreement with previous findings [47,48], showing no bias for heterotetramer formation in CHO cells expressing a mixture of WT and G219P mutant NaChBac monomers and dimers. hyperpolarization-activated whole-cell currents [47], (Supplementary Figure S1A,B). The mutation shifts the voltage of the maximal current from −10 mV in WT to −180 mV in the mutant (Supplementary Figure S1D), thus currents at −10 mV originating from separate channel populations of WT and L226P homotetramers can be easily separated. The current recordings from CHO cells cotransfected with NaChBac-encoding WT:L226P cDNAs in a ratio 3:1 (Supplementary Figure S1C,E) exhibited unique currents at −10 mV, which cannot be explained by the simple addition of whole current traces from homotetramer channels formed from either L226P or wild-type NaChBac (note that there is no current at −10 mV from L226P channels), indicating that unique heterotetramers are being formed. Assuming that the assembly of heterotetramers is formed without bias, the proportions of channel types can be determined by binomial distribution. It is noteworthy that this assumption is in agreement with previous findings [47,48], showing no bias for heterotetramer formation in CHO cells expressing a mixture of WT and G219P mutant NaChBac monomers and dimers.

*Entropy* **2020**, *22*, x FOR PEER REVIEW 6 of 20

channel gating of NaChBac from depolarization-activated whole-cell currents to non-inactivating

Using this approach, CHO cells were co-transfected with cDNAs of NaChBac-encoding WT and mutants, with varied Q*<sup>f</sup>* in the SF, in different ratios. Note that the open probabilities and single channel conductances for the WT NaChBac (LESWAS) and LEDWAS homotetramers were equivalent (Supplementary Figure S2), and that the whole-cell Na<sup>+</sup> currents from cells expressing homotetramer WT and LEDWAS were similar in magnitude (Figure 1a,c), consistent with the expression of the channel (i.e., number of channels) being independent of the single amino acid mutations introduced into the SF. Figure 1 shows the currents recorded from cells transfected with defined mixtures of NaChBac monomers; see Table 1 for the probabilities of different charged species, assuming that the assembly follows a binomial distribution. Using this approach, CHO cells were co-transfected with cDNAs of NaChBac-encoding WT and mutants, with varied Q*f* in the SF, in different ratios. Note that the open probabilities and single channel conductances for the WT NaChBac (LESWAS) and LEDWAS homotetramers were equivalent (Supplementary Figure S2), and that the whole-cell Na+ currents from cells expressing homotetramer WT and LEDWAS were similar in magnitude (Figure 1a,c), consistent with the expression of the channel (i.e., number of channels) being independent of the single amino acid mutations introduced into the SF. Figure 1 shows the currents recorded from cells transfected with defined mixtures of NaChBac monomers; see Table 1 for the probabilities of different charged species, assuming that the assembly follows a binomial distribution.

**Figure 1.** Na+/Ca2+ selectivity for NaChBac monomer mixtures. The voltage–current relations (**A**,**B**) and the mean (+/− SEM) whole-cell peak current density at −10 mV (**C**,**D**) for Na+ (**A**,**C**) and Ca2+ (**B**,**D**) in CHO cells transfected with cDNAs encoding for NaChBac channels possessing either a wild-type selectivity filter (LESWAS/LES) or a mutated selectivity filter (LASWAS/LAS or LEDWAS/LED); 5 μg of total DNA was used per transfection and was composed of either a mixture of types of cDNA at **Figure 1.** Na+/Ca2<sup>+</sup> selectivity for NaChBac monomer mixtures. The voltage–current relations (**A**,**B**) and the mean (+/<sup>−</sup> SEM) whole-cell peak current density at <sup>−</sup>10 mV (**C**,**D**) for Na<sup>+</sup> (**A**,**C**) and Ca2<sup>+</sup> (**B**,**D**) in CHO cells transfected with cDNAs encoding for NaChBac channels possessing either a wild-type selectivity filter (LESWAS/LES) or a mutated selectivity filter (LASWAS/LAS or LEDWAS/LED); 5 µg of total DNA was used per transfection and was composed of either a mixture of types of cDNA at defined ratios, as indicated in Table 1 and on the *X*-axis, or a single cDNA type. Numbers in parentheses indicate the number of replicates.



Whole-cell currents were initially recorded in bath solution containing 140 mM of Na-methanesulfonate, followed by recordings after the complete replacement of bath Na<sup>+</sup> with 100 mM of Ca-methanesulfonate. Na<sup>+</sup> permeation appears relatively insensitive to changes in Q*<sup>f</sup>* values between −4*e* and −8*e* and equivalent in cells expressing only LESWAS and/or LEDWAS monomers. Focusing on channels exhibiting Q*<sup>f</sup>* values less than −4*e*, it is interesting to note that, despite co-transfection with a 1:1 ratio of LASWAS:LESWAS resulting in an expected only 6.2% of the channel population being homotetramers of LESWAS (Q*<sup>f</sup>* = <sup>−</sup>4*e*), the Na<sup>+</sup> current density was approximately 30% of that recorded from cells expressing only LESWAS homotetramer channels (Figure 1a,c). An equivalent interpretation can be made for measurements of current density from cells transfected with a 1:3 ratio of cDNAs encoding LASWAS:LESWAS: the sodium current density was equivalent to that recorded from cells expressing only LESWAS homotetramers, despite only 32% of the channel population being predicted to be homotetrameric LESWAS. It is also interesting to compare the current density of the 3:1 LASWAS:LESWAS-expressing cells. Note that these cells show about 25% current density compared to the LESWAS-only cells (5 and 20 pA/pF, respectively). If one looks at the binominal predictions, 25% of channels are predicted to have Q*<sup>f</sup>* = −2*e* and greater and this is consistent with a Q*<sup>f</sup>* = −1*e* and 0 being non-conducting (Figure 1c). The simplest explanation for the disproportionately large Na<sup>+</sup> current in cells expressing mixtures of LESWAS and LASWAS monomers is that functional NaChBac channels possessing a SF with a Q*<sup>f</sup>* value less than −4*e* are functional and able to mediate the Na<sup>+</sup> influx.

Extending this type of analysis to the Ca2<sup>+</sup> currents, cells transfected with a 1:3 ratio of LESWAS:LEDWAS encoding cDNAs (in which 0.3% of expressed functional channels were predicted to be LEDWAS homotetramers, with a Q*<sup>f</sup>* <sup>=</sup> <sup>−</sup>8*e*) exhibited a similar current density for Ca2<sup>+</sup> influx as that from cells expressing only LEDWAS channels (Figure 1d). Thus, functional NaChBac channels possessing SFs with a Q*<sup>f</sup>* value of less than <sup>−</sup>8*<sup>e</sup>* appear to be able to mediate Ca2<sup>+</sup> influx, with the possibility that a Q*<sup>f</sup>* value of <sup>−</sup>5*<sup>e</sup>* is sufficient to permit Ca2<sup>+</sup> permeation. This explanation is also consistent with the observation that the Ca2<sup>+</sup> current density is greatest in cells transfected with equal and 1:3 ratios of LESWAS:LEDWAS (Figure 1b,d). Note that the Na<sup>+</sup> influx current density remains relatively constant in cells transfected with both LESWAS and LEDWAS encoding cDNAs, indicating that the effect of varying Q*<sup>f</sup>* between <sup>−</sup>4*<sup>e</sup>* and <sup>−</sup>8*<sup>e</sup>* was specific to the Ca2<sup>+</sup> current density.

### 3.1.2. Na+/Ca2<sup>+</sup> Selectivity for Concatenated NavMS Channels

Although the use of a mixed population of cDNAs encoding for NaChBac and its mutants suggested the value of Q*<sup>f</sup>* to be a major determining factor for Na+/Ca2<sup>+</sup> selectivity, the results are subject to the caveat that the whole-cell currents result from the cumulative current from an unknown but predictable range of different channel types. To address this complication, we attempted to generate a stable concatenation of NaChBac to enable the expression of a homogeneous population of NaChBac mutants; however, we have previously shown [38] the NaChBac oligomer to be unstable and not to remain intact in the plasma membrane. In contrast, an equivalent intact NavMs oligomer could be stably expressed in HEK293T cells [38] and thus enable the generation of a homogeneous population of bacterial channels, in which the Q*<sup>f</sup>* value of the SF can be altered in steps of 1*e*. The SF of eukaryotic CaVs is formed by a ring of glutamates (the EEEE locus) and a conserved aspartate residue in domain II (D2p51 [37]). The D2p51 residue is suggested to form a binding site for a third incoming Ca2<sup>+</sup> from the extracellular side of the pore and thus bring an additional positive charge to the SF region necessary for the release of a bound Ca2<sup>+</sup> to the cytosolic side (i.e., a knock-on mechanism [49]). Although direct evidence for the role of the D2p51 in Ca2<sup>+</sup> permeation remains elusive, replacing the D2p51 residue in Cav1.2 (aka D707) with a neutral amino acid residue significantly reduces the Ca2<sup>+</sup> binding of the SF [29]. Thus, to gain further insight into the role of the D2p51 in Ca2<sup>+</sup> permeation, we used site-directed mutagenesis targeted to repeat I or II in the NavMs oligomer to generate a bacterial NaV with an "EEEED" locus (Q*<sup>f</sup>* = −5*e*) in the SF (Supplementary Figure S3). NavMs has a

high homology (45% sequence identity) to NaChBac [28,45], which should enable comparison with the results from NaChBac. *Entropy* **2020**, *22*, x FOR PEER REVIEW 10 of 20

The WT NavMs SF is defined by <sup>177</sup>LESWSM182, and we generated NavMs tetramers (Supplementary Figure S3) with the S179D mutation in either repeat I (mutant DI) or repeat II (mutant DII). Both mutants are therefore expected to carry a charge of −5*e* in the SF. Figure 2 shows typical whole-cell currents from WT and mutant NavMs in bath solution containing either 140 mM of Na<sup>+</sup> or 100 mM of Ca2+. The WT NavMs SF is defined by 177LESWSM182, and we generated NavMs tetramers (Supplementary Figure S3) with the S179D mutation in either repeat I (mutant DI) or repeat II (mutant DII). Both mutants are therefore expected to carry a charge of −5*e* in the SF. Figure 2 shows typical whole-cell currents from WT and mutant NavMs in bath solution containing either 140 mM of Na+ or 100 mM of Ca2+.

**Figure 2.** Na+/Ca2+ selectivity for NavMs concatemer possessing varied Q*f* values in their SF. The original recordings representatives of wild-type NavMS (**A**) and its DI (**B**) and DII (**C**) mutants in 140 mM Na+ solution (grey traces) and in 100 mM Ca2+ solution (black traces). The voltage–current relations (**D**,**E**) and the mean (+/− SEM) whole-cell peak current density at −10 mV (**F,G**) for Na+ (**D,F**) and Ca2+ (**E,G**) in HEK 293T cells transfected with cDNAs encoding for either wild-type or mutated NavMS. Numbers in parentheses indicate the number of replicates; and in HEK293T cells transfected with wild-type NavMS concatemer (WT) and mutant NavMS concatemer (DI and DII). **Figure 2.** Na+/Ca2<sup>+</sup> selectivity for NavMs concatemer possessing varied Q*<sup>f</sup>* values in their SF. The original recordings representatives of wild-type NavMS (**A**) and its DI (**B**) and DII (**C**) mutants in 140 mM Na<sup>+</sup> solution (grey traces) and in 100 mM Ca2<sup>+</sup> solution (black traces). The voltage–current relations (**D**,**E**) and the mean (+/<sup>−</sup> SEM) whole-cell peak current density at <sup>−</sup>10 mV (**F,G**) for Na<sup>+</sup> (**D,F**) and Ca2<sup>+</sup> (**E,G**) in HEK 293T cells transfected with cDNAs encoding for either wild-type or mutated NavMS. Numbers in parentheses indicate the number of replicates; and in HEK293T cells transfected with wild-type NavMS concatemer (WT) and mutant NavMS concatemer (DI and DII).

In order to make quantitative comparisons between the electrophysiological behavior of NaChBac and NavMs mutants, the peak calcium and sodium currents as well as their ratio are tabulated in Supplementary Table S3 for the NaChBac heterotetramer populations and in Supplementary Table S2 for the NavMs mutants. The comparison of the data of the two tables shows that the ratio of the peak current densities for Na+ and Ca2+ in wild-type NavMs (0.018) is comparable to that for wild-type NaChBac (0.010). The tables also show that the ratio of peak current densities for Ca2+ and Na+ in the two NavMs mutants with an SF charge of −5*e* (0.080 for the DI mutant and 0.097 for the DII mutant) is similar to that for mutant channels formed from the expression of the 3LES:1LED mixture of NaChBac (0.054) in CHO cells, which yields the highest probability of occurrence (42%) of heterotetramers with an SF charge equal to −5*e*. Although both data sets support the increased Ca2+ permeability in −5*e* mutant bacterial sodium channels, the difference in the Ca2+ In order to make quantitative comparisons between the electrophysiological behavior of NaChBac and NavMs mutants, the peak calcium and sodium currents as well as their ratio are tabulated in Supplementary Table S3 for the NaChBac heterotetramer populations and in Supplementary Table S2 for the NavMs mutants. The comparison of the data of the two tables shows that the ratio of the peak current densities for Na<sup>+</sup> and Ca2<sup>+</sup> in wild-type NavMs (0.018) is comparable to that for wild-type NaChBac (0.010). The tables also show that the ratio of peak current densities for Ca2<sup>+</sup> and Na<sup>+</sup> in the two NavMs mutants with an SF charge of −5*e* (0.080 for the DI mutant and 0.097 for the DII mutant) is similar to that for mutant channels formed from the expression of the 3LES:1LED mixture of NaChBac (0.054) in CHO cells, which yields the highest probability of occurrence (42%) of heterotetramers with an SF charge equal to <sup>−</sup>5*e*. Although both data sets support the increased Ca2<sup>+</sup> permeability in <sup>−</sup>5*<sup>e</sup>* mutant bacterial sodium channels, the difference in the Ca2<sup>+</sup> current magnitude that is evident upon

that factors other than the value of Q*f* are important in determining the Ca2+ permeability.

current magnitude that is evident upon comparing NavMs and NaChBac channels clearly indicates

comparing NavMs and NaChBac channels clearly indicates that factors other than the value of Q*<sup>f</sup>* are important in determining the Ca2<sup>+</sup> permeability.

#### *3.2. Computational Results*

In an attempt to gain a molecular-level understanding of the different behavior of WT NavMs (Q*<sup>f</sup>* = −4*e*) and its mutant with charge Q*<sup>f</sup>* = −5*e*, we ran equilibrium MD simulations in a 100 mM solution of CaCl<sup>2</sup> or 140 mM NaCl (for 150 and 100 ns, respectively). The initial structure of WT NavMs was taken from the Protein Data Bank (ID: 3ZJZ). Mutation S179D on chain A and embedding in a membrane of 248 POPC molecules was performed using the CHARMM membrane builder [40,41].

In 140 mM of NaCl solution, the WT NavMs SF is stably occupied by a single Na<sup>+</sup> even if transient events of occupation by a second ion can also be spotted (Supplementary Figure S4a). Q*<sup>f</sup>* = −5*e* mutants SF is almost immediately occupied by two Na<sup>+</sup> ions and, after 40 ns, the filter becomes stably occupied by three sodium ions (Supplementary Figure S4b). The different behavior of the two species is also reflected in the PMF profile (Figure 3a,b), which is characterized by a single deep minimum centered at z = 4–5 Å for WT NavMs, and a minimum split into three sub-basins at z = 2.0 Å, z = 5.0 Å, and z = 8–9 Å, corresponding to three different binding regions, for the Q*<sup>f</sup>* = −5*e* mutant. The barriers between the sub-basins are in the order of 1–2 kcal/mol and can be easily overcome at the simulation temperature, yet the sodium ions linger in each binding site for longer than they would in case of a uniform probability distribution of occupancy.

The nature of these binding sites can be better characterized by analyzing the conformation of the SF in the last frame of the 100 ns simulations (Figure 3e–h). A notable feature of wild-type NavMs is that the side chains of E178 residues do not point toward the center of the channel, but they are aligned along the channel wall pointing towards the extracellular side. As a result, the distance between the resident sodium ion and the ε-oxygen of E178 always exceeds 4.0 Å. This means that there are no direct sodium-protein interactions; Na<sup>+</sup> interacts with the protein via water molecules in its hydration shell. Indeed, the withdrawn placement of E178 side chains leaves sufficient space in the SF for Na<sup>+</sup> to fully keep its first hydration shell of six water molecules. In contrast, the conformation of the SF of Q*<sup>f</sup>* = −5*e* mutant revealed three sodium ions that directly interact with the residues of the SF; specifically, the extracellular one interacting with D179 and E178 both located on chain A, the central one with E178, and the intracellular sodium with the backbone carbonyl group of one of the L177 residues. The additional negative charge thus determines an enhanced ability of the NavMs mutant to capture sodium ions from the bulk. This, combined with the possibility of a knock-on mechanism deriving from the simultaneous presence of three Na<sup>+</sup> ions in the SF, possibly explains the larger sodium current density for NavMs channels with Q*<sup>f</sup>* = −5*e*. As a result of this structural arrangement (and in contrast to that for the WT; Figure 3c), Na<sup>+</sup> ions accessing the SF of the mutant lose on average 3.5 water molecules. However, this loss is compensated by the interactions with the oxygens provided by the acidic residues (2 oxygens) and by other protein residues (1 oxygen), such that the total number of coordinating oxygens is maintained equivalent to that for sodium in bulk solution (Figure 3c,d). Note that the interactions of the resident ions with all other residues of the SF are water-mediated.

In 100 mM of CaCl<sup>2</sup> within the timescale of our simulations, no Ca2<sup>+</sup> gains access to the SF of the WT channel, while a single Ca2<sup>+</sup> enters into the SF of the Q*<sup>f</sup>* <sup>=</sup> <sup>−</sup>5*<sup>e</sup>* mutant during the early stages of the simulation and thereafter remains locked inside, while also repelling other potentially incoming calcium ions (Supplementary Figure S5a,b). This pattern is in keeping not only with the behavior of WT and mutant NavMs concatemer, but also with the results of the experiments on mixed populations of NaChBac heterotetramers. In fact, while the calcium peak current of the LESWAS homotetramers (Q*<sup>f</sup>* = −4*e*) is just 0.66 pA/pF, that of the 3LES:1LED population, where we expect the highest proportion of channels with Q*<sup>f</sup>* = −5*e*, is tenfold higher (6.0 pA/pF). The seeming mismatch between the currents recorded in experiments and the total block of the Ca2<sup>+</sup> ion revealed by the simulations is (at least in part) due to the fact that, in the latter, no electric field was applied. Moreover, experimental recordings are performed on a timescale of hundreds of milliseconds, one million times longer than that covered

by simulations, allowing time for slow, activated events of ion permeation. A comparison between our computational results and those by other groups is discussed in Supplementary Figure S6. *Entropy* **2020**, *22*, x FOR PEER REVIEW 12 of 20

**Figure 3.** MD simulations for the WT and mutant NavMS in NaCl 140 mM. (**a**,**b**) Potential of Mean Force of Na+ as a function of the axial position in WT NavMS (**a**) and the mutant with charge Q*f* = −5*e* (**b**). (**c**,**d**) Average number of coordinating oxygens per sodium ion in axial bins with a thickness of 2.0 Å. (**c**) Wild-type NavMs; (**d**) NavMs mutant with Q*f* = −5*e*. The distance cutoff to identify sodiumchloride interactions was set to 3.5 Å, and for sodium-oxygen to 3.2 Å. Color code is as follows. Blue line: number of coordinating water-provided oxygens; green line: number of coordinating oxygens provided by aspartate and glutamates; red line; number of coordinating oxygens provided by other protein residues; black line: total number of coordinating oxygens; magenta line: number of coordinating chlorides. (e,h) Configuration of the selectivity filter of wild-type NavMS (**e**,**g**) and the mutant with charge Q*f* = −5*e* (**f**,**h**). All the structures correspond to the last frame of a 100 ns simulation in 0.14 M NaCl. Panels (**e**,**f**) show a side view of the SF; panels (**g**,**h**) show the top view. Glu178 is shown in red, while Asp179 is shown in purple. The backbone of Leu177 is shown in green. Sodium ions are portrayed as blue beads. Panels (**e**,**g**) also show the water molecules that mediate the interactions between the resident sodium ion and the protein in wild-type NavMS. The nature of these binding sites can be better characterized by analyzing the conformation of **Figure 3.** MD simulations for the WT and mutant NavMS in NaCl 140 mM. (**a**,**b**) Potential of Mean Force of Na<sup>+</sup> as a function of the axial position in WT NavMS (**a**) and the mutant with charge Q*<sup>f</sup>* <sup>=</sup> <sup>−</sup>5*<sup>e</sup>* (**b**). (**c**,**d**) Average number of coordinating oxygens per sodium ion in axial bins with a thickness of 2.0 Å. (**c**) Wild-type NavMs; (**d**) NavMs mutant with Q*<sup>f</sup>* = −5*e*. The distance cutoff to identify sodium-chloride interactions was set to 3.5 Å, and for sodium-oxygen to 3.2 Å. Color code is as follows. Blue line: number of coordinating water-provided oxygens; green line: number of coordinating oxygens provided by aspartate and glutamates; red line; number of coordinating oxygens provided by other protein residues; black line: total number of coordinating oxygens; magenta line: number of coordinating chlorides. (**e**,**h**) Configuration of the selectivity filter of wild-type NavMS (**e**,**g**) and the mutant with charge Q*<sup>f</sup>* = −5*e* (**f**,**h**). All the structures correspond to the last frame of a 100 ns simulation in 0.14 M NaCl. Panels (**e**,**f**) show a side view of the SF; panels (**g**,**h**) show the top view. Glu178 is shown in red, while Asp179 is shown in purple. The backbone of Leu177 is shown in green. Sodium ions are portrayed as blue beads. Panels (**e**,**g**) also show the water molecules that mediate the interactions between the resident sodium ion and the protein in wild-type NavMS.

the SF in the last frame of the 100 ns simulations (Figure 3e–h). A notable feature of wild-type NavMs

The position of the ion in the SF revealed by the Potential of Mean Force (PMF) shows that Ca2<sup>+</sup> ions do visit the vestibule region of the WT channel, but they never enter into the SF (Figure 4a). The presence of an additional negative charge in the SF (S179D) is sufficient to pull in a Ca2<sup>+</sup> ion that occupies a binding site centered at z = 6.5 Å (Figure 4b). The PMF minimum corresponding to this binding site has a depth of approximately 7.0 kcal/mol, which, at the simulation temperature of 300 K, corresponds to 11.5 *kBT*. The energy well is thus so deep that a single Ca2<sup>+</sup> cannot leave the SF. Thus, to be consistent with the experimental recording of Ca2<sup>+</sup> current (Figure 2e,g), calcium permeation must involve some sort of knock-on mechanism. The role of the aspartate residue in the SF is immediately highlighted by Figure 4e–h, which shows the configuration of the SF in the last frame of the simulation. The resident calcium ion appears to be directly bonded to D179 and to E178, both located on chain A (Figure 4). The interactions with the other glutamates of the SF are all water-mediated. In order to better characterize calcium hydration, in Figure 4c,d we plot the average number of coordinating oxygen atoms per calcium ion in axial bins with a thickness of 2.0 Å. Figure 4 shows that when a calcium ion enters into the SF, the number of hydrating water molecules drops from approximately 8 to 4.5. This dehydration is compensated by an increase in the number of coordinating oxygens provided by aspartate and glutamate residues (approximately 3). Thus, when Ca2<sup>+</sup> enters the SF, the total number of coordinating oxygens remains roughly unchanged (Figure 4c,d).

A collective diffusion model approach was adopted to approximate the Ca2<sup>+</sup> currents [46]. The algorithm relates the spontaneous permeation events at equilibrium with steady currents induced by small voltages. This approach thus enables the estimation of currents from equilibrium simulations; however, as it is based on linear response theory, its predictions are reliable only in a small voltage range. The results of the calculation are summarized in Table 2.

**Table 2.** Current estimates through linear response theory. The first column shows the NavMs species analyzed, either the wild-type form EEEE with an SF charge Q*<sup>f</sup>* = −4*e* or the mutant EEEED with an additional negative charge in the SF (Q*<sup>f</sup>* = −5*e*). The second column shows the ion carrying the current, the third column reports the estimated conductance in pS, and the fourth column lists the estimated current at V = −20 mV. This voltage corresponds to the peak current in the current–voltage plots determined from whole-cell patch-clamp experiments.


Notwithstanding the limitations of our calculations, the collective diffusion modelling predictions are in reasonable agreement with the experimental observations (Figure 2). For example, (1) whole-cell recordings showed that peak sodium currents increased by approximately 1.5-fold in the Q*<sup>f</sup>* = −5*e* NavMs channel (−35 to −55 pA/pF); this is mirrored by a 1.5-fold increase in sodium conductance (from 23.06 to 35.37 pS) predicted by linear response theory calculations. (2) Experimental measurements of peak calcium currents in the Q*<sup>f</sup>* = −5*e* mutant are approximately 10 times smaller than those for sodium. This is consistent with the modelling in the Q*<sup>f</sup>* = −5*e* mutant, in which a seven-fold greater sodium (35.37 pS) conductance is predicted compared to that for calcium (4.87 pS). (3) The small finite Ca2<sup>+</sup> influx predicted in the WT NavMs (Table 2; 1.69 pS) can be observed in the electrophysiological recordings (Figure 2e).

*Entropy* **2020**, *22*, x FOR PEER REVIEW 14 of 20

**Figure 4.** MD simulations for WT and mutant NavMs in CaCl2 100 mM. (**a**,**b**) Potential of Mean force of Ca2+ as a function of the axial position in WT NavMS (**a**) and the mutant with charge Q*f* = −5*e* (**b**). (**c**,**d**) Average number of coordinating oxygens per calcium ion in axial bins with a thickness of 2.0 Å. (**c**) Wild-type NavMs; (**d**) NavMs mutant with Q*f* = −5*e*. The distance cutoff to identify both calciumchloride and calcium-oxygen interactions was set to 3.5 Å. Color code is as follows. Blue line: number of coordinating water oxygens; green line: number of coordinating oxygens provided by aspartate and glutamates; red line; number of coordinating oxygens provided by other protein residues; black line: total number of coordinating oxygens; magenta line: number of coordinating chlorides. (**e**–**h**) Configuration of the selectivity filter of wild-type NavMS (**e**,**g**) and the mutant with charge Q*f* = −5*e* (**f**,**h**). All the structures correspond to the last frame of a 150 ns simulation in 0.10 M of CaCl2. Panels (**e**,**f**) show a side view of the SF; panels (**g**,**h**) show the top view. Glu178 is shown in red, while Asp179 is shown in purple. The backbone of Leu177 is shown in green. Calcium ions are portrayed as orange beads. **Figure 4.** MD simulations for WT and mutant NavMs in CaCl<sup>2</sup> 100 mM. (**a**,**b**) Potential of Mean force of Ca2<sup>+</sup> as a function of the axial position in WT NavMS (**a**) and the mutant with charge Q*<sup>f</sup>* <sup>=</sup> <sup>−</sup>5*<sup>e</sup>* (**b**). (**c**,**d**) Average number of coordinating oxygens per calcium ion in axial bins with a thickness of 2.0 Å. (**c**) Wild-type NavMs; (**d**) NavMs mutant with Q*<sup>f</sup>* = −5*e*. The distance cutoff to identify both calcium-chloride and calcium-oxygen interactions was set to 3.5 Å. Color code is as follows. Blue line: number of coordinating water oxygens; green line: number of coordinating oxygens provided by aspartate and glutamates; red line; number of coordinating oxygens provided by other protein residues; black line: total number of coordinating oxygens; magenta line: number of coordinating chlorides. (**e**–**h**) Configuration of the selectivity filter of wild-type NavMS (**e**,**g**) and the mutant with charge Q*<sup>f</sup>* = −5*e* (**f**,**h**). All the structures correspond to the last frame of a 150 ns simulation in 0.10 M of CaCl<sup>2</sup> . Panels (**e**,**f**) show a side view of the SF; panels (**g**,**h**) show the top view. Glu178 is shown in red, while Asp179 is shown in purple. The backbone of Leu177 is shown in green. Calcium ions are portrayed as orange beads.

A collective diffusion model approach was adopted to approximate the Ca2+ currents [46]. The algorithm relates the spontaneous permeation events at equilibrium with steady currents induced by small voltages. This approach thus enables the estimation of currents from equilibrium simulations;

#### **4. Discussion and Conclusions**

In this work, we engineered radial asymmetry in the bacterial NaChBac and NavMs channels as a first attempt to mimic the features of eukaryotic voltage-gated sodium and calcium channels. It is well known that prokaryotic sodium channels are characterized by a glutamate ring that imparts a charge <sup>−</sup>4*e* to the SF and endows the channel with Na<sup>+</sup> selectivity. It is also well established that an increase in the negative charge of the SF makes the channel progressively more calcium-selective. Pioneering studies by the Clapham group, for instance, showed that mutating into aspartate either serine of the SF sequence TLESWAS of NaChBac decreases the PNa/PCa ratio, while a mutation of both serines makes the channel completely calcium-selective [20]. Using the same strategy, more recently Tang et al. replaced the TLESWSM sequence in the SF of NavAb with TLDDWSD, causing a complete shift from sodium to calcium selectivity [21]. It is noteworthy that, due to the tetrameric symmetry of prokaryotic NaVs, in all these studies the charge of the SF was always varied in −4*e* steps and radial symmetry was maintained. It is thus known that a charge Q*<sup>f</sup>* = −4*e* is typical of a Na+-selective channel, while a charge −8*e* or −12*e* leads to calcium selectivity. This change in the Q*<sup>f</sup>* value is rather coarse and does not address the fact that the SF of eukaryotic channels is asymmetric. Therefore, it is important to investigate the influence on the selectivity of charge changes by −1*e* steps.

The study of random heteroteramers in our work indicated that channels with an SF charge smaller than <sup>−</sup>4*e* mediate Na<sup>+</sup> currents and channels with an SF charge in the <sup>−</sup>4*e* < Q*<sup>f</sup>* < <sup>−</sup>8e range conduct Ca2+. Furthermore, the study of the NavMs concatemer showed that the presence of an additional negative charge in the SF leads to a significant increase in the Na<sup>+</sup> and Ca2<sup>+</sup> current.

The electrophysiological behavior of the NaChBac populations of randomly assembled heterotetramers appears to be in reasonable agreement with the predictions of the Ionic Coulomb Blockade (ICB) model [23–25]. According to this model, ion permeation and selectivity through channels mainly depend on the Q*<sup>f</sup>* of the SF. If calcium permeation is plotted as a function of Q*<sup>f</sup>* , a pattern of alternating conductance and stop bands can be observed. In contrast, the same plot for sodium predicts a steep increase in current magnitude up to values in Q*<sup>f</sup>* of <−2*e*, followed by a plateau and the absence of stop bands (Figures 2 and 3, in [25]), consistent with sodium permeation being relatively insensitive to changes in Q*<sup>f</sup>* . Thus, the predictions of the ICB model appear to be compatible with the plot of peak sodium currents in Figure 1c. Furthermore, it is tempting to envisage the pattern of Ca2<sup>+</sup> current density shown in Figure 1d as an oscillation in the calcium conductance (i.e., conductance and stop bands), which would repeat over a wider range of Q*<sup>f</sup>* values. The ICB model, however, appears to be less successful in explaining the so-called "EEEE paradox"—that is, the apparently shared "EEEE" motif in both the sodium-selective bacterial NaVs and the calcium-selective eukaryotic CaVs. Kaufman et al. tentatively reconciled this inconsistency, noting the presence of a conserved aspartate close to the EEEE ring of eukaryotic CaVs, thus redefining the motif as EEEED, which raises the SF charge to −5*e* [36]. Our experiments on the NavMs concatemer go some way towards confirming this prediction, but also highlight that other factors in addition to the value of Q*<sup>f</sup>* are important. The ICB model predicts that a charge <sup>−</sup>5*<sup>e</sup>* allows the access of a third Ca2<sup>+</sup> ion when the SF is already occupied by two resident calcium ions. Our MD simulations, however, show that while no calcium ion gains access to the SF of WT NavMs, only a single Ca2<sup>+</sup> ion stably occupies the filter of the mutant with charge <sup>−</sup>5*e*. This calcium ion is strongly bound to D179 and E178 located on the same subunit, and sits in a free energy well so deep that it cannot leave the SF. At the same time, the resident ion probably exerts an electrostatic repulsion on other potentially incoming Ca2<sup>+</sup> ions, preventing a knock-on mechanism in a similar fashion as that described for NaChBac [49].

Our experiments and simulations thus suggest that the extra negative charge is effective in the capture of cations from the bulk, but it does not promote permeation. Contrary to that postulated by simplified physical models (in which the channel atomic structure is not considered), the charge of the SF is not the only determinant of conduction and selectivity. It is possible that calcium flow in eukaryotic Cavs requires some sort of fine modulation of charge effects. Flood et al., for instance, performed an interesting computational study grafting the SF and external vestibular region of the

human NaV1.2 channel into the scaffold of the NavRh bacterial channel [35]. Their multi-microsecond MD simulations revealed that permeation and selectivity depend on the close interplay of the DEKA and EEDD rings, so that the charge of the extended filter region is −5*e* as in our NavMs mutant. In its protonated state, the lysine residue of the DEKA ring acts like a built-in sodium ion involved in the formation of multi-carboxylates/multi-ion complexes. When the charged ammonium group of lysine is in the HFS site, where the electrostatic potential is most negative, it creates a smooth electrostatic environment leading into the cavity, whereas, when it is bent toward the central cavity, it creates a zone of high electrostatic potential that cuts the cavity off from the SF. Our recent work [38], showing the possibility to create stable concatemers of the bacterial NavMs channel, offers the opportunity to experimentally test these computational predictions by creating a bacterial channel chimera where the SF and vestibule of the human Nav1.2 channel are grafted onto the NavMs concatemer.

Since no positively charged residue appears to be located close to the SF of eukaryotic CaVs, the fine modulation of the charge might rely on the differential protonation of the acidic residues of the EEEED locus. The effect of protonation has been extensively studied through MD simulations. Furini et al., for instance, showed that the glutamate side chains in NavAb can adopt two different orientations pointing either towards the extracellular environment or towards the central cavity [34]. Interestingly, they found that the likelihood of the inwardly directed arrangement increases when E177 residues are protonated. Moreover, the presence of a glutamate residue with the side chain directed to the central cavity increases the energy barrier for the translocation of sodium ions. Since E177 was observed to adopt an alternative conformation in MD simulations with Ca2<sup>+</sup> ions [50], it is possible that these protonation-induced configurations also affect selectivity. While the control of the protonation state of the filter is a trivial task in MD simulations, it is a challenging endeavor in biophysical experiments.

This leads us to the methodological aspect of our work. Our study not only tested the importance of SF charge in controlling ion selectivity and permeation, but created new tools extending the use of bacterial channels as models of eukaryotic ones. Indeed, the current work is the first one to report experiments on a NaV channel in which the pore region has been mutated to have radial asymmetry, and thus it represents an important first step in bridging the major limitation in using bacterial sodium channels to investigate their eukaryotic counterparts. Our methodology will enable us to design physical experiments to investigate the mechanisms of the fine modulation of charge effects that are likely to occur in asymmetric eukaryotic channels, such as that predicted by Flood et al. [35].

A further methodological merit of our approach is its relevance in understanding the effect of pH on channel permeation and selectivity. In fact, when the pH is varied, the four glutamates of the SF are unlikely to be protonated or deprotonated simultaneously. A more probable scenario is that they are protonated or deprotonated one at a time, resulting in +1*e* or −1*e* changes in the SF charge [23]. Finally, our combination of molecular dynamics and electrophysiological approaches provided fresh insight into the molecular mechanisms of cation permeation in bacterial sodium channels, and gave insight into understanding the molecular mechanisms that underlie the function of NaVs and CaVs.

**Supplementary Materials:** The following are available online at http://www.mdpi.com/1099-4300/22/12/1390/s1: Supplementary Table S1 Parameters of NavMs equilibration. Supplementary Table S2. Peak currents ratios of WT and mutant NavMs channels. Supplementary Table S3. Peak currents ratios of NaChBac heterotetramers. Supplementary Figure S1. The tetramer formation via co-transection is proved to be a random process without bias for homo- or hetero-tetramer formation. Supplementary Figure S2. The single channel currents recorded from WT NaChBac. Supplementary Figure S3. Schematic representation of NavMs concatemer. Supplementary Figure S4. Ion occupancy of Selectivity Filter in MD simulations of WT NavMs and its mutant with Q*<sup>f</sup>* = −5*e* in NaCl 140 mM. Supplementary Figure S5. Ion occupancy in Selectivity Filter in MD simulations of WT NavMs and its mutant with Q*<sup>f</sup>* = −5*e* in CaCl<sup>2</sup> 100 mM. Supplementary Figure S6. Current-voltage plot calculations: comparison of constant electric field simulations and equilibrium simulations in conjunction with linear response theory.

**Data Availability:** Data related to this research are openly available from the University of Warwick archive at (https://wrap.warwick.ac.uk/143573). Fedorenko, Olena A., Khovanov, Igor A., Roberts, Stephen K., and Guardiani, Carlo (2020) Data for *Changes in ion selectivity following asymmetrical addition of charge to the selectivity filter of bacterial sodium channels* [Dataset].

**Author Contributions:** Conceptualization, O.A.F., C.G., S.K.R., and I.A.K.; methodology, O.A.F. and C.G.; formal analysis, O.A.F. and C.G.; investigation, O.A.F. and C.G.; writing—original draft preparation, O.A.F. and C.G.; writing—review and editing, S.K.R. and I.A.K. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by EPSRC; project: Ionic Coulomb blockade oscillations and the physical origins of permeation, selectivity, and their mutation transformations in biological ion channels (grant numbers: EP/M015831/1 and EP/M016889/1). CG is currently supported by a project that has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No. 803213).

**Acknowledgments:** We are grateful to Huaping Sun for the generation of the −5*e* NavMs constructs.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

#### *Perspective*

## **Prospects of Observing Ionic Coulomb Blockade in Artificial Ion Confinements**

#### **Andrey Chernev, Sanjin Marion and Aleksandra Radenovic \***

Laboratory of Nanoscale Biology, Institute of Bioengineering, School of Engineering, EPFL, 1015 Lausanne, Switzerland; andrey.chernev@epfl.ch (A.C.); sanjin.marion@epfl.ch (S.M.)

**\*** Correspondence: aleksandra.radenovic@epfl.ch

Received: 9 November 2020; Accepted: 11 December 2020; Published: 18 December 2020

**Abstract:** Nanofluidics encompasses a wide range of advanced approaches to study charge and mass transport at the nanoscale. Modern technologies allow us to develop and improve artificial nanofluidic platforms that confine ions in a way similar to single-ion channels in living cells. Therefore, nanofluidic platforms show great potential to act as a test field for theoretical models. This review aims to highlight ionic Coulomb blockade (ICB)—an effect that is proposed to be the key player of ion channel selectivity, which is based upon electrostatic exclusion limiting ion transport. Thus, in this perspective, we focus on the most promising approaches that have been reported on the subject. We consider ion confinements of various dimensionalities and highlight the most recent advancements in the field. Furthermore, we concentrate on the most critical obstacles associated with these studies and suggest possible solutions to advance the field further.

**Keywords:** nanofluidics; ionic Coulomb blockade; 2D materials; nanopores; nanotubes; angstrom slits

#### **1. Introduction**

In the past fifteen years, various artificial nanofluidic platforms have become highly compelling for fundamental studies of physical phenomena and numerous practical applications where the transport of the confined ions plays a crucial role [1,2]. Among the most exciting practical applications are power generation [3–7], filtration, and molecular separation [8–10]. The last five years have witnessed remarkable progress in the fabrication of nanofluidic devices, enabling researchers to develop artificial nanofluidic systems with the confinement of one to a few water molecules (below 1 nm) [11]. Such nanofluidic platforms have been realized in zero-dimensional (0D), 1D, or 2D geometry [7,12–14] (Figure 1**)**. These platforms exhibit giant permeability and ion selectivity comparable to biological ion channels, excluding anions and macromolecules, and closely mimic functionalities previously observed only in biological channels. Superb selectivity of the sodium/calcium family of channels has fascinated scientists due to its physiological relevance and underlying physical mechanism [15–17]. A simplified electrostatic and Brownian dynamics model of the prototypical model Ca2<sup>+</sup> or Na<sup>+</sup> channel has been used in describing its conduction and selectivity [18], echoing the phenomenology of Coulomb blockade [17]. The term ionic Coulomb blockade (ICB) was suggested first as a counterpart of the electronic Coulomb blockade (ECB) by Krems and Di Ventra [19]. Confinement below 1 nm dictates a departure from the mean-field assumptions as the correlations between the ions and finite-size effects cannot be neglected and complicates the insight into the ionic charge transfer at the nanoscale [17,20,21].

*Entropy* **2018**, *20*, x 2 of 14

**Figure 1.** Different geometries are considered for confined nanofluidic systems from two‐dimensional (2D) nanoslits to 1D nanowires and quantum dot‐like nanopores in the atomically thin material matrix, i.e., 2D nanopores. **Figure 1.** Different geometries are considered for confined nanofluidic systems from two-dimensional (2D) nanoslits to 1D nanowires and quantum dot-like nanopores in the atomically thin material matrix, i.e., 2D nanopores.

Artificial nanofluidic devices developed for single‐ion conductivity measurements provide a unique opportunity to test the proposed models and reveal the many‐body effects in ionic systems. This perspective article aims to discuss practical challenges in verifying the models that evoke ionic Coulomb blockade in a variety of settings that are now available due to recent progress in nanofabrication [12,22,23]. Furthermore, it also provides suggestions for the integration of additional approaches such as modulation of the charges, pressure, pH, ionic strength, temperature, and potential in the range not attainable on lipid bilayers. Artificial nanofluidic devices developed for single-ion conductivity measurements provide a unique opportunity to test the proposed models and reveal the many-body effects in ionic systems. This perspective article aims to discuss practical challenges in verifying the models that evoke ionic Coulomb blockade in a variety of settings that are now available due to recent progress in nanofabrication [12,22,23]. Furthermore, it also provides suggestions for the integration of additional approaches such as modulation of the charges, pressure, pH, ionic strength, temperature, and potential in the range not attainable on lipid bilayers.

As suggested initially, ionic Coulomb blockade is based on the relation of ion kinetic and barrier energies and manifests as the nonlinear transport of ions. This behavior was later experimentally observed in atomically thin, sub‐nanometer‐sized molybdenum disulfide (MoS2) nanopores [13]. Although theoretical models predict ICB and suggest particular conditions to reveal the effect [17,20,24], irrefutable experimental observation of this phenomenon is still challenging. The reason for that is a whole set of issues that have to be taken into account before considering this complex process. The most critical issues that could mask ICB observation are related to the fact that ionic transport is measured at room temperatures, causing large charge fluctuations, which leads to the increased noise, and instability of the nanofluidic devices. Together with regular wetting and contamination challenges in nanofluidics [25], we can argue that the ionic systems' nonlinear current‐ voltage characteristics are insufficient to prove the ICB effect. Unlike in the case of ECB [26], there are still no convincing data showing conductance oscillations or single‐ion devices. As suggested initially, ionic Coulomb blockade is based on the relation of ion kinetic and barrier energies and manifests as the nonlinear transport of ions. This behavior was later experimentally observed in atomically thin, sub-nanometer-sized molybdenum disulfide (MoS2) nanopores [13]. Although theoretical models predict ICB and suggest particular conditions to reveal the effect [17,20,24], irrefutable experimental observation of this phenomenon is still challenging. The reason for that is a whole set of issues that have to be taken into account before considering this complex process. The most critical issues that could mask ICB observation are related to the fact that ionic transport is measured at room temperatures, causing large charge fluctuations, which leads to the increased noise, and instability of the nanofluidic devices. Together with regular wetting and contamination challenges in nanofluidics [25], we can argue that the ionic systems' nonlinear current-voltage characteristics are insufficient to prove the ICB effect. Unlike in the case of ECB [26], there are still no convincing data showing conductance oscillations or single-ion devices.

Therefore, it is essential to critically assess the suitability of different nanofluidic platforms for the observation of ICB and propose novel solutions that will allow us to untangle the role of the issues mentioned above. New geometries like nanopores in 2D materials, single‐digit nanotubes and angstrom slits play a key role and are essential in verifying the different parameters required to observe ICB reproducibly. The complex nature of ionic transport can be considered with different settings and thus allow direct observation of experimental parameter variation such as pore length, diameter, access resistance, surface charge, gating voltage temperature, and many others [27]. Therefore, it is essential to critically assess the suitability of different nanofluidic platforms for the observation of ICB and propose novel solutions that will allow us to untangle the role of the issues mentioned above. New geometries like nanopores in 2D materials, single-digit nanotubes and angstrom slits play a key role and are essential in verifying the different parameters required to observe ICB reproducibly. The complex nature of ionic transport can be considered with different settings and thus allow direct observation of experimental parameter variation such as pore length, diameter, access resistance, surface charge, gating voltage temperature, and many others [27].

#### **2. Prerequisites for Ionic Coulomb Blockade**

The basic prerequisites for ionic Coulomb blockade are shared with its electronic counterpart [28,29]. The geometrical requirement for ICB to occur is represented by two reservoirs with charge carriers, in our case ions, with a channel in between with a large resistance to ion transport at least in one direction [17,19,20,24]. More specifically, the channel has to provide the type of confinement that allows ions to dwell inside of it as it was suggested originally [15]. Ion transport between the two reservoirs is inhibited due to strong Coulomb interactions and an ion dwelling at or near the channel, so that there is a limit on the possible charge, which can be transferred between the two sides of the chamber [17,19,20,24]. When the energy barrier for ions to enter the channel, and subsequently traverse it, is larger than the system's thermal energy, ionic Coulomb blockade can be observed. We identify two main mechanisms that are expected to result in ionic Coulomb blockade signatures in nanopore or nanochannel systems. One is linked to capacitive charging [19,30] and analogous to the electronic case of a tunnel junction. The second mechanism is linked to the existence of an "island" corresponding to a quantum dot in ECB, which can be gated [15,16,20] and can be further developed towards the analogy of the single-electron transistor [20].

The first case, which we call the capacitive ICB, is relatively similar to the ECB case with a single tunneling junction. It occurs if an ion transits between the two reservoirs and charges it analogously to a capacitor, thus producing a barrier for further ions to transit [19,30].

This case requires that the channel/pore be ion-selective so that a capacitive barrier can form for an ion of valence *z*, that the thermal energy is lower than the capacitive self-energy of the channel *U* = *Q*2/2*C<sup>s</sup>* , and that the transferred charge *Q* = *nze* dwells inside the pore or next to the pore entrance for a sufficiently long time for it to interact with the transport of further ions. We note that the ICB effect is expected to be stronger for ions carrying more than one unit of charge, as the self-energy has a quadratic dependence on the valence *z.* The capacity of the channel/pore of length *L* and radius *R* to store charge is the self-capacitance *C<sup>s</sup>* and can, neglecting any fine effects, be approximated by *C<sup>s</sup>* = 4π0*rR* <sup>2</sup>/*L*, with *<sup>r</sup>* the relative dielectric constant of water. In general, the capacitive self-energy is a good measure if ICB is detectable as it gives information on the resulting energy barriers for ion transport.

Coarse-grained molecular dynamics simulations [19,30] indicate capacitive ICB could occur in sub-nm pores in quasi-2D or 2D materials but also predict a weak nonlinear current-to-voltage dependence uncharacteristic for the classical electronic Coulomb blockade. In contrast, Brownian dynamics simulations show that the role of screening and ion-pair formation via Bjerrum pairs are paramount in these systems [20]. Here, both parameters are strongly modulated by the dielectric constant of water, which has been shown to be reduced in confinement in respect to its bulk value [31,32], and which would decrease the capacitive self-energy and require lower temperatures or smaller and thinner pores for ICB to be detectable. Moreover, finite-size effects would need to be considered, including the peculiarities of pore-wall interactions, which can strongly influence ion transport in quasi-2D or 2D membranes [33–36].

In the second case, ionic Coulomb blockade happens in the systems that are reminiscent of a quantum dot, confined with two energy barriers—schemes that are used for ECB observations and applications. As for the ICB counterpart, it can be observed in the setting where the ion occupancy of the channel limits ion transport [17,20] and where it has a role analogous to an isolated quantum dot in the electronic case [29]. Here, ions need to have an energetically favorable position inside the channel, which binds them there in spite of thermal motion so that their presence would block other ions traversing the channel. This stable position can be formed due to ion-channel interactions, either by electrostatic gating voltage modulation or by the presence of surface charges [24,37]. In the latter case, the surface charge will then attract oppositely charged ions and thus block ion transport through the channel.

When a fixed voltage bias is applied between the reservoirs, the ion current oscillates with peaks at certain quantized values of these gating surface charges. A typical signature of ICB is when this neutralizing charge has a stepwise dependence on the gating voltage or surface charge, and when the conductance of the channel has peaks at certain values with conductance almost completely suppressed otherwise [17].

Microscopically, there are several distinctions between the ionic and the electronic Coulomb blockades. While electrons are negatively charged by default, ions can be both negative and positive, with different valences available. These ions are also present in a solution so that the electrostatic interaction between them is screened, losing its long-range components. Thus, such an effect will be local and will require extreme confinements for size exclusion and Coulomb interactions to dominate transport properties. In the case of strong 1D confinement, ion transport through such channels is expected to go through the dissociation of pairs of ions (second Wien effect) [38,39].

Nevertheless, in the case of gating via voltage or wall charge inhibiting the movement of ions due to tight binding, some of the ions are tightly bound, and only under certain conditions are they able to dissociate, a situation recently termed as a "fractional Wien effect" [20]. Under this model, it is also predicted that when electrostatic screening becomes too large, the ionic Coulomb blockade effect becomes suppressed. This also provides an upper limit on the size of a confining channel to about 2.5 nm. This means that, in certain geometries, ICB could be detected even in pores larger than one nanometer, depending on the magnitude of electrostatic screening.

Providing unquestionable proof for the ionic Coulomb blockade is not an easy task. Most experimental "wet" implementations cannot easily change the interaction between the ions and the channel interior, i.e., by changing its surface charge or applying a gating voltage to it. One of the major indicators of ICB is nonlinear current-to-voltage curves [17,20,24,35]. As the bias voltage is increased until a certain threshold is crossed, there is no ion transport through the system. After that, a nonlinear increase in current is expected until linear (ohmic) behavior is obtained. In the classical electronic CB, several steps called the Coulomb staircase in the current versus voltage curves are expected in this nonlinear regime [29]. Recent work indicates that, in the case of ICB, there could be no such steps, with a direct transition to an ohmic regime [20]. However, it is problematic to identify ICB only with current-voltage characteristics (IV curves), especially in the absence of a staircase-like pattern. Realistic devices are known to exhibit leakage currents that will exceed the current in the blocked ICB conditions [13], possibly overshadowing any conduction steps.

Furthermore, similar nonlinear IV curves are known to occur during the electrowetting of channels and are caused by remnant gas inside the pore system [25,40–42]. Partial wetting of these systems presents with nonlinear activation-like IV curves, similar to the ones that have been associated with ionic current blockade [13]. There are many open questions regarding our understanding of ICB, which could help to design better experiments. Variations of the ambient temperature could be used to probe the activation energy barriers and possibly promote ICB. Varying the solution's pH could be used to change the gating charges, thus enabling the variation of the gating charge in larger channels. However, the main requirement is still to manufacture a device with such properties to optimize the conditions for the observation of ICB, a topic we tackle in the following sections.

#### **3. ICB in 2D Nanopores**

Over the last decade, 2D materials have become a rich area of research and are showing tentative signs of impacting our everyday life [43]. In bulk, these materials have the form of layered crystals, with van der Waals interaction holding together 2D layers with a thickness starting from 0.3 nm. These 2D sheets are under intense study because of their fascinating electronic properties, spanning the range from isolating and semiconducting to superconducting. Nanopores formed in 2D membranes from graphene [44–46], hexagonal boron nitride (h-BN) [47], transition metal dichalcogenides (TMDCs) [4,48], and MXenes [49] have been used to investigate nanofluidic phenomena or can act as a single molecular sensor.

2D nanopores have the peculiarity that due to their thickness, or lack of it, the pore sizes need to be small for the interactions between a single ion and the pore to dominate ion transport. For capacitive ICB to occur, we need ion selectivity, a sufficient capacitive barrier, and ions to remain next to the pore entrance for a sufficient time. 2D nanopores are known for their selectivity to ions, which is expected to come from the surface charges and the resulting ion enhancement zones near the pore [5]. Simple classical estimates of the self-capacitance would support 2D material nanopores having a sufficient energy barrier [13], but there is no evidence of ion retention for sufficiently long times to cause this blockade effect unless the ion is retained inside the pore itself. Classical bulk models for nanopore conductance break down as we approach the nanometer scale [21,33,50], requiring us to include the specifics of the pore atomic structure and its interactions with the ions and water.

Ion transport in 2D systems is a rich topic of research [33,51], but little is known experimentally about the peculiarities of transport in sub-nm pores. When the degree of confinement approaches the ion's size, i.e., sub-nm size ranges, then an effect called hydration layer shredding occurs [19,37,52]. If an ion needs to lose a part of its hydration layer to traverse the pore, it experiences an energy barrier that can also result in a nonlinear IV curve, one of the expected signatures for ICB. Unless there is some sort of weakly bound state for the ion in the middle of the channel, this hydration layer shredding is not enough by itself to induce a Coulomb blockade effect. This can be achieved by tailoring pore interior bond edges, which could temporarily bind an ion inside the pore, causing it to block further ion transport [50]. Pore edge termination is known to influence ion transport [53–55] and could be tailored by binding other chemical species and modifying the pore interior [54–58]. Ion transport through pores is determined by bulk solution and solvent-mediated interactions between the ions and the pore interior edges. This interaction has been shown to cause ions to dwell inside the pore [24,35]; however, it is unclear if this influences ion transport via the Coulomb blockade mechanism.

2D nanopore systems have several unresolved issues that still require further study. Control of the pore shapes and sizes is not easy to achieve as pores are known to be highly sensitive to the fabrication method [47,59] (Figure 2). However, this can be addressed by the predefined pore shapes. Together with electrostatic gating, it is expected to be a highly effective approach for ion selectivity and controlled transport [24]. In this case, the trapped cations themselves not only act as natural clogs of the pores, but also "protect" the clogged pore states by creating the repulsive potential around each pore that essentially suppresses the knock-on mechanism by mobile cations dancing around.

Solvation of pore systems is expected to be able to change the pore structure, and etching of pores and their enlargement due to chemical species or dissolved reactive oxygen. In addition, 2D systems are prone to surface contaminants (e.g., hydrocarbons) that originate from the 2D material transfer process [60]. These contaminants can clog pores or change their properties, making controlled experiments difficult. Additionally, wetting or clogging issues can plague these systems [25] and possibly even produce nonlinear IV curves due to the presence of nanobubbles [40,41,61]. Moreover, the standard method of wetting all nanoscale devices uses prewetting with water-ethanol mixtures, a procedure known to produce nanobubbles on hydrophobic surfaces and which has been linked to the same effect on mildly hydrophobic 2D material surfaces [25]. In addition, remnants of bound ethanol, or other molecules, have been shown to modify surface properties and can persist after subsequent flushing [62]. Differentiating these effects from ICB nonlinear IV curves requires additional probes and separate gating control.

*Entropy* **2018**, *20*, x 6 of 14

**Figure 2.** (**a**–**c**) Selection of the most popular 2D materials, i.e., graphene, hexagonal boron nitride, and molybdenum disulfide structures in bulk and monolayer form and generation and evolution of defects in 2D materials under e‐beam. Defects and pores are created under 80 keV e‐beam irradiation while imaging with a transmission electron microscope. Images adapted from (**a**) [59] and (**b**,**c**) [48]. **Figure 2.** (**a**–**c**) Selection of the most popular 2D materials, i.e., graphene, hexagonal boron nitride, and molybdenum disulfide structures in bulk and monolayer form and generation and evolution of defects in 2D materials under e-beam. Defects and pores are created under 80 keV e-beam irradiation while imaging with a transmission electron microscope. Images adapted from (**a**) [59] and (**b**,**c**) [48].

Solvation of pore systems is expected to be able to change the pore structure, and etching of pores and their enlargement due to chemical species or dissolved reactive oxygen. In addition, 2D systems are prone to surface contaminants (e.g., hydrocarbons) that originate from the 2D material transfer process [60]. These contaminants can clog pores or change their properties, making controlled experiments difficult. Additionally, wetting or clogging issues can plague these systems [25] and possibly even produce nonlinear IV curves due to the presence of nanobubbles [40,41,61]. Moreover, the standard method of wetting all nanoscale devices uses prewetting with water‐ethanol mixtures, a procedure known to produce nanobubbles on hydrophobic surfaces and which has been linked to the same effect on mildly hydrophobic 2D material surfaces [25]. In addition, remnants of bound ethanol, or other molecules, have been shown to modify surface properties and can persist after subsequent flushing [62]. Differentiating these effects from ICB nonlinear IV curves requires additional probes and separate gating control. Gating can be applied to nanofluidic platforms for modulating ICB in 2D pores (Figure 3). Gating can be applied to nanofluidic platforms for modulating ICB in 2D pores (Figure 3). Notably, it can be applied in various ways, and not only limited to electrostatic gating. For example, one approach would be applying mechanical stress directly to nanopores [63] as the strain has been predicted to modulate ion transport through sub-nm 2D nanopores by promoting or removing a stable dwelling spot for ions in the center of the pore [31–33]. Strain moves the pore edge atoms out of their equilibrium positions, changing ion–atom and ion-water interactions, thus modifying the free energy landscape for the ions, possibly promoting or destroying a bound state for the ion inside the pore. Chemical gating of these systems is limited, and the surface charge can be modulated via pH modulation [11] and photo-gating [6], but it is not expected to make significant changes to the interaction of ions with the pore edges and could promote pore growth due to etching. Another useful probe would be studying the effect of temperature variation to provide information about the energy barriers for transport and possibly promote or suppress ICB effects. Although promising, this strategy has a reduced temperature range, in essence, from the melting to the boiling point of the solution.

Notably, it can be applied in various ways, and not only limited to electrostatic gating. For example, one approach would be applying mechanical stress directly to nanopores [63] as the strain has been predicted to modulate ion transport through sub‐nm 2D nanopores by promoting or removing a stable dwelling spot for ions in the center of the pore [31–33]. Strain moves the pore edge atoms out of their equilibrium positions, changing ion–atom and ion‐water interactions, thus modifying the free

pore. Chemical gating of these systems is limited, and the surface charge can be modulated via pH modulation [11] and photo‐gating [6], but it is not expected to make significant changes to the interaction of ions with the pore edges and could promote pore growth due to etching. Another useful probe would be studying the effect of temperature variation to provide information about the energy barriers for transport and possibly promote or suppress ICB effects. Although promising, this

**Figure 3.** Common 2D nanopore measurement system with temperature control (**a**) equipped with transverse contacts for additional measurement channel or gating (**b**) and 2D material encapsulation (hBN‐MoS2‐hBN) (**c**) for the improved signal‐to‐noise ratio and decoupling of the measurement **Figure 3.** Common 2D nanopore measurement system with temperature control (**a**) equipped with transverse contacts for additional measurement channel or gating (**b**) and 2D material encapsulation (hBN-MoS<sup>2</sup> -hBN) (**c**) for the improved signal-to-noise ratio and decoupling of the measurement circuits (cis-trans and transverse).

#### circuits (cis‐trans and transverse). *Heterostructures*

*Heterostructures* A promising approach may be to use layered heterostructures of 2D materials [64]. Here, a stable position for an ion inside the pore could be designed through a combination of different materials in A‐B or A‐B‐A scheme where the monolayers are held together by van der Waals forces (Figure 3c). The major challenge for this is the existence of hydrocarbon contaminants complicating the reproducible production of such heterostructures using material transfer, where one 2D material (B) is transferred over another (A). This disadvantage would not be an issue in the case of in situ growth of heterostructures. Similar heterostructures have already been opened in several different research avenues, such as graphene interacting with h‐BN allowed several groups to study the Hofstadter butterfly effect while numerous optoelectronic devices were based on the heterostructures from A promising approach may be to use layered heterostructures of 2D materials [64]. Here, a stable position for an ion inside the pore could be designed through a combination of different materials in A-B or A-B-A scheme where the monolayers are held together by van der Waals forces (Figure 3c). The major challenge for this is the existence of hydrocarbon contaminants complicating the reproducible production of such heterostructures using material transfer, where one 2D material (B) is transferred over another (A). This disadvantage would not be an issue in the case of in situ growth of heterostructures. Similar heterostructures have already been opened in several different research avenues, such as graphene interacting with h-BN allowed several groups to study the Hofstadter butterfly effect while numerous optoelectronic devices were based on the heterostructures from semiconducting monolayers [65–67].

semiconducting monolayers [65–67]. Van der Waals heterostructures could be used to achieve direct gating by sandwiching a semiconductor such as MoS2 in between two layers of a wider bandgap energy 2D material such as h‐BN. This suggested geometry could be considered either on a stand‐alone basis or in combination with the in‐plane transverse gating similar to the previous nanopore field‐effect transistor (FET) Van der Waals heterostructures could be used to achieve direct gating by sandwiching a semiconductor such as MoS<sup>2</sup> in between two layers of a wider bandgap energy 2D material such as h-BN. This suggested geometry could be considered either on a stand-alone basis or in combination with the in-plane transverse gating similar to the previous nanopore field-effect transistor (FET) approaches [68,69].

approaches [68,69]. Another promising setting could be achieved by using stacked membranes (i.e., graphene oxide), where the percolating path of ions through the porous membranes and between the layers Another promising setting could be achieved by using stacked membranes (i.e., graphene oxide), where the percolating path of ions through the porous membranes and between the layers could be modified due to ions getting stuck and blocking ion transport [70].

could be modified due to ions getting stuck and blocking ion transport [70]. However, to progress in the aforementioned systems, one has to have solid fabrication approaches for each of them. It is important to advance in nanopore fabrication techniques to reliably However, to progress in the aforementioned systems, one has to have solid fabrication approaches for each of them. It is important to advance in nanopore fabrication techniques to reliably produce sub-nm pores through layered 2D membranes using existing techniques, namely an electron beam in transmission electron microscope [50,60,61,71] or electron beam lithography [72,73]. This approach is known to be the most common practice in nanopore research. However, complex multiple-step processes often result in relatively high contamination and low device yield; thus, alternative solutions have to be applied.

One possible solution is a 2D material growth over the aperture in a suspended SiN membrane [73]. However, the first attempts show that the SiN pore ends up inducing growth, and thus, the island

of growth is located exactly above it. Another fine point is that the resulting 2D material might be polycrystalline, which makes its structure controversial for pore stability. Following that, this approach requires a significant improvement in order to benefit from the potential cleanliness of such a process.

On the contrary, the electrochemical reaction allows the precise control of the pore size in the flow cell [74]. However, this technique does not allow the pore size confirmation other than using ionic current-voltage characteristics. The reason for that is a poor quality of nanopore samples after the nanofluidic experiments, as they appear to be irreversibly contaminated after the measurements in a flow cell.

#### **4. ICB in 1D Nanowires**

The 1D channels represent the next promising nanofluidic platform (Figure 4). Such systems, reminiscent of narrow nanotubes, have been proposed as a very effective ion-confining solution capable of showing the single-ion transport and consequently able to reveal an ICB [20]. Furthermore, based on the geometrical confinement properties, Fermi-Dirac distribution for ions can be achieved, and thus, single-ion transport was predicted [17] and explained in detail for biological channels [17] and infinite *Entropy* 1D channels [ **2018**, *20*, 20x ]. 9 of 14

**Figure 4.** 1D experimental geometry suggested for the ionic Coulomb blockade observation (**a**–**c**) in 1D nanochannels and (**d**) 2D nanoslits by the application of an external electric field and a consequent electrostatic gating of the ion‐confining channel. Pictures (a–c) reprinted from Kavokine et al. [20]. **Figure 4.** 1D experimental geometry suggested for the ionic Coulomb blockade observation (**a**–**c**) in 1D nanochannels and (**d**) 2D nanoslits by the application of an external electric field and a consequent electrostatic gating of the ion-confining channel. Pictures (**a**–**c**) reprinted from Kavokine et al. [20].

One possible solution to that would be an extensive application of the most advanced fabrication techniques allowing the electrostatic control over each channel **(**Figure 4d**)** and thus, limiting the ionic current to a single channel of interest, a similar level of control could also be achieved mechanically similarly as in microfluidics [79]. Moreover, the 2D nanoslit fabrication approach [12] lacks control over channel side wall quality, and despite the atomically smooth surface of top and bottom of the channel, this might significantly affect the transport in narrow nanoslits with the comparable aspect ratio in both horizontal and vertical directions. Nonlinear transport of ions in such systems has been shown [14]; however, in that case, it is not associated with the Coulomb blockade of the channels due to a multiple input from a selection of nanoscale effects that might affect the conductance [80]. One possible solution to that would be an improved electrostatic control over the electric field in the channel as mentioned earlier (Figure 4d), similarly as proposed for other geometries. **6. Summary** State-of-the-art techniques that allow single nanotube manipulation [3,75], if applied to the ultra-narrow nanotubes, could allow fabrication of a single conductive channel for ions, 1D ionic channels that show great potential for the observation of ICB. Inspired by ECB applications [76], 1D ionic channels were also suggested as the charge-carrier pumping systems, where the gating voltage on multiple gating contacts oscillates out of phase and thus allows a turnstile mechanism and a single-ion passage through the channel [20]. The unifying point of the aforementioned theoretical works [17,20] is the fixed wall charge that has to be controlled in order to manipulate the charge carrier transport at a single-ion level. However, despite the great success in nanofabrication [3,75], theoretically predicted effects are still to be confirmed experimentally, partly due to an even higher complexity of nanofabrication required (point contact, low leakage current, room temperature of measurements, high device yield). However, from the technical point of view, nanofabrication approaches developed for carbon nanotubes and nanowire field-effect transistors could be of help when applied to one-dimensional ion channels [75,77].

#### The complex nature of ion transport in nano‐confined systems is usually challenged with natural artifacts, such as improper wetting of the pore, pore instability (especially in the case of a 2D **5. ICB in 2D Nanoslits**

to be fulfilled, i.e., *kT<<Ec.*

the pore on a single‐ion level.

confining geometry (either electrical or optical).

nanopore), inability to confirm the pore size after the measurement, and many others. Development of a well‐defined approach that would be able to address these issues one by one would propel the First of all, 2D nanoslits inherited the flexibility of material stacking that was suggested for the first time in graphene research [66]. Following that, as in the case with 2D nanopores, this type of ion

Temperature control is important in such systems since it has a significant impact on the ions' mobility, the potential barrier size, and charge fluctuations, which might increase noise and the instability of nanofluidic devices. Therefore, in both ICB and ECB cases, the temperature criterion has

Electrostatic or optical gating of the confinement is expected to provide a crucial transition from a closed to an opened state within the system. The current status of nanofabrication allows us to create wafer‐scale batches of devices showing a corresponding quality [22], and to combat one of the most important bottlenecks in nanopore research field—the small amount of nanopore devices and consequently low statistics of the measurements. Together with the demonstrated pioneer approaches [68,69,81] the nanopore research field is reaching the limitations and heading towards the achievement of the controlled confinement that would allow us to probe charged states within

However, we assume that it is important to show the most controversial points that still have to be resolved in the field. For example, nanoscale contamination is starting to play an even more crucial role in nanoelectronic devices. Leakage currents and electrochemical reactions happening at the nanoscale in 2D materials and gating contacts of nanopore field‐effect transistors come into play

research field. However, this conception sounds much easier than it is; one would think of the

confinement in 2D nanoslits could benefit from a wide variety of available materials [59]. Furthermore, 2D nanoslits attract a particular interest as a first-ever system to experimentally show the value of the dielectric constant of water in confinement [32], frictionless water transport [78], and unmatched selectivity, allowing only water and protons to pass through the narrow 2D nanoslits [7]. However, despite the significant flexibility given by the fabrication process of 2D nanoslits, one particularly important property has to be considered—these systems are represented by multiple nanofluidic channels, and thus, read-out of a single-channel ion transport and control over the transport at a single-ion level is still complicated.

One possible solution to that would be an extensive application of the most advanced fabrication techniques allowing the electrostatic control over each channel (Figure 4d) and thus, limiting the ionic current to a single channel of interest, a similar level of control could also be achieved mechanically similarly as in microfluidics [79]. Moreover, the 2D nanoslit fabrication approach [12] lacks control over channel side wall quality, and despite the atomically smooth surface of top and bottom of the channel, this might significantly affect the transport in narrow nanoslits with the comparable aspect ratio in both horizontal and vertical directions. Nonlinear transport of ions in such systems has been shown [14]; however, in that case, it is not associated with the Coulomb blockade of the channels due to a multiple input from a selection of nanoscale effects that might affect the conductance [80]. One possible solution to that would be an improved electrostatic control over the electric field in the channel as mentioned earlier (Figure 4d), similarly as proposed for other geometries.

#### **6. Summary**

The complex nature of ion transport in nano-confined systems is usually challenged with natural artifacts, such as improper wetting of the pore, pore instability (especially in the case of a 2D nanopore), inability to confirm the pore size after the measurement, and many others. Development of a well-defined approach that would be able to address these issues one by one would propel the research field. However, this conception sounds much easier than it is; one would think of the experiments where the local temperature is controlled, together with the further gating of the confining geometry (either electrical or optical).

Temperature control is important in such systems since it has a significant impact on the ions' mobility, the potential barrier size, and charge fluctuations, which might increase noise and the instability of nanofluidic devices. Therefore, in both ICB and ECB cases, the temperature criterion has to be fulfilled, i.e., *kT*<<*Ec.*

Electrostatic or optical gating of the confinement is expected to provide a crucial transition from a closed to an opened state within the system. The current status of nanofabrication allows us to create wafer-scale batches of devices showing a corresponding quality [22], and to combat one of the most important bottlenecks in nanopore research field—the small amount of nanopore devices and consequently low statistics of the measurements. Together with the demonstrated pioneer approaches [68,69,81] the nanopore research field is reaching the limitations and heading towards the achievement of the controlled confinement that would allow us to probe charged states within the pore on a single-ion level.

However, we assume that it is important to show the most controversial points that still have to be resolved in the field. For example, nanoscale contamination is starting to play an even more crucial role in nanoelectronic devices. Leakage currents and electrochemical reactions happening at the nanoscale in 2D materials and gating contacts of nanopore field-effect transistors come into play [68,69]. With irreproducible point-like defects and contaminants, ICB observation is even more challenging since, unlike in ECB, one cannot simply "freeze" the unwanted states and free charges. This problem keeps growing as long as we need more fabrication steps to produce a better-controlled ion-confining system.

Furthermore, recent studies [25] have shown that nanoscale contamination may lead to partial rewetting, which adds a high level of complexity for these studies. In particular, we have to admit that the most common practice of pore wetting with ETOH solution with subsequent electrolyte flushing and exchange may lead directly to the formation of nanoscale gas bubbles [61]. Apart from that, the typical ETOH wetting process has been shown to change the hydrogen distribution over the surface, and thus, change the ion-solvation part of the interactions of the ions with the surface and consequently mask the ICB.

Nevertheless, we would like to draw attention to the most promising future prospects and applications of this field to probe each contribution to ion transport phenomena and to improve the control over these platforms. That would allow us to achieve a synthetic model of a cell's ion channels on a discrete channel level and benefit from step-by-step additions to a solid-state platform's controls. As a result, we may expect to obtain crucial insights on gating, strain, temperature control, etc., to use as inputs in the bottom-up computation design of biopores [82], which have already proven themselves as a unique platform for sensing and single-molecule experiments [83]. So far, no further ICB signatures have been experimentally demonstrated as conductance oscillations and Coulomb staircase, even though these have been predicted for biological nanopores [17].

**Author Contributions:** Writing-Original Draft Preparation, A.C. and S.M.; Writing-Review & Editing, A.C., S.M., A.R.; Visualization, A.R.; Supervision, A.R.; Project Administration, A.R.; Funding Acquisition, A.R. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by Swiss National Science Foundation (SNSF) support through grant number [200021\_192037]" and "The A.C. was funded by Eurotech Postdoc Programme (European Commission Horizon 2020. Grant Agreement number 754462).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Entropy* Editorial Office E-mail: entropy@mdpi.com www.mdpi.com/journal/entropy

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18

www.mdpi.com ISBN 978-3-0365-1645-5