Brownian Motion Paving the Way for Molecular Translocation in Nanopores

,


Introduction
[8][9] Electroosmosis in an SSNP is primarily a consequence of the presence of surface charge on the dielectric nanopore-sidewall; [10] the surface charge induces the formation of an electrical double layer (EDL) immediately above the surface of the nanoporesidewall hosting counterions of the surface charge. [11]Application of an external transmembrane electric voltage causes these counterions in the EDL to drift along the electric field. [12,13]The movement of the counterions sets the nearby water molecules in motion due to drag force (because of viscosity) thereby generating the electroosmotic flow (i.e., electroosmosis).Electroosmosis is more pronounced in smaller nanopores wherein the EDL constitutes a larger fraction of their smallest constriction; EDL is ≈1 nm in thickness for physiologically relevant electrolytes. [10]eometrically asymmetric SSNPs have been shown to rectify not only ionic current [14,15] but also molecular translocation. [8,9]ectification refers specifically to higher currents or more frequent translocations in one direction than in the other of an SSNP.By clarifying the causal chain connecting the key physical factors and processes leading to rectification, an analytical model has recently been proposed to account for the rectified ionic current in asymmetric SSNPs as well as in their sibling nanopipettes. [16]19][20][21] However, the analysis is generally challenged by uncertainties in capturing the rapidly translocating analytes [22] or registering distorted signals [8] that may also be mistaken as noise.By assuming an amplifier with a bandwidth of 1 MHz and then a temporal resolution of 1 μs, only 6% of the total protein translocation events could be missed according to a first-passage time distribution model. [23]Thus, it is necessary to develop higher-bandwidth measurement systems to proceed with detection regardless of the type of proteins.A 250 kHz measurement system corresponding to the minimum dwell time of ≈2.5 μs has been reported to distinguish proteins with small molecular weights below 30 kDa. [22] Moreover, a solid-state sensing platform with an extraordinary bandwidth of ≈10 MHz has been shown to be able to register a very short dwell time of ≈0.2 μs in single-strained DNA translocation events. [24]High bandwidth with the advantage of capturing signals from fast-passing analytes of diminishing dwell time is, unfortunately, accompanied by rapidly increasing noise with frequency in the high-frequency end above 10 kHz. [25,26]Hence, the actual detectable dwell time can be compromised.
Instead of striving for higher bandwidth electronics, studying the rectification of analyte translocations could be an alternative to exploring translocation behavior and nanopore design to advance analyte sensing based on SSNPs.The rectification could be represented by the FTE ratio defined as the quotient of the measured mean FTE at negative bias to that at positive bias. [8,9]The translocation data of streptavidin (molecular weight ≈60 kDa) and IgG 1 (≈150 kDa), using a measurement platform with 10 kHz bandwidth, show a clear dependence of FTE ratio on the protein size relative to the smallest constriction of a truncated-pyramidal nanopore (TPP). [8]This dependency and the observation of an unprecedentedly strong rectification of DNA translocation in a bowl-shaped nanopore (BNP) [9] have been analyzed based on numerical simulation of electroosmosis.Unfortunately, this numerical approach only yields indirect information inferred from ion flow patterns and risks missing critical details pertaining to the analyte translocation itself especially considering the limited knowledge about analyte movements due to the poor temporal resolution of the measurement platform used.
The present work investigates the translocation dynamics based on COMSOL Multiphysics under conditions mimicking the actual experiment configurations as closely as possible.Studying the trajectory of each and every analyte in a large analyte cohort assuming the basic properties of streptavidin remains a focus.It is well known that for analytes of 1 μm or below in diameter, the influence of Brownian forces on their moving trajectories becomes significant.The influence is especially pronounced for analytes of a few nanometers in diameter with nanoporetranslocation times in μs or shorter using SSNPs.A key development of the model implemented on COMSOL is the consideration of Brownian force in addition to the already established electroosmotic and electrophoretic forces.Motivated by previous analytical studies based on, e.g., continuous random Gaussian field, [27] we have developed a numerical simulation platform by including Brownian motion of the analytes in the presence of electrical field and fluid.The trajectories of a large analyte cohort translocating a nanopore are simultaneously monitored to yield statistically meaningful translocation behaviors.The simulation confirms some previously unknown, yet critical, subtleties with respect to analyte translocation in nanopores.The Brownian motion was previously invoked as an alternative to dielectrophoresis for particle-particle interaction at a distance in micro-/nanofluidics. [28] The simulation in terms of the FTE ratio compares qualitatively well with the experiment, despite the significant difference in time scale between the simulation (10 ms) and experiment (100 s) as well as the simplifications necessary for the model development and implementation.Details revealed by the simulation are implicative for the design of new SSNP sensors.

Geometrical Shape
The foremost factor for performing the simulation is to set up an appropriate geometry along with material properties in the model to be implemented.In terms of geometry, the shape of truncated conical nanopores (TCPs) was adopted with similar dimensions to those of the experimental TPP.The structural details of such a TCP are illustrated in Figure 1a-d to best mimic the experimental setup around TPP.To enable resource-demanding computation, the 2D axis-symmetric system is divided into upper and lower reservoirs relative to the conical nanopore located in the middle of the system.Each reservoir has a diameter of 2 μm and a height of 1 μm (Figure 1c).Simultaneously, the morphology of the nanopore was set to be the same as in the experiment (details can be found in Tables S1 and S2 Supporting Information).

Physical Properties
The nanopore region was set to Si with a native oxide layer (≈1.5 nm) and the rest to water (Figure 1d).The COMSOL material library was utilized to retrieve the representative parameters of water, i.e., density, viscosity, and relative permittivity (Table S1 Supporting Information).Models to account for the dependence of viscosity on size [29] and salt concentration [30,31] were also incorporated.A recent study shows that the relative viscosity of nanoconfined water tends to increase with a decrease in the size of nanotubes. [29]The resultant relative viscosity of ≈2.3 to ≈1.1 depends on the height (z) and hence on the horizontal radius of the TCP (Note S1 and Figure S1 Supporting Information).As done previously, [8] the surface charge density of the naturally oxidized Si nanopore sidewall was set to −0.02 C m −2 to approximate the experimentally determined result of −0.016 C m −2 .Besides, induced surface charge (ISC) due to the external electric field, [32] whose density is position-dependent along the nanopore sidewall, was also considered in the model.

Ionic Movement
The simulation of analyte translocations was divided into two steps: to calculate electroosmotic flow in the fluid region and to track the time evolution of input spheres.The ionic movement was jointly described by the Poisson equation, the Nernst-Planck equation, and the Navier-Stokes equation, where, Φ, F, , z i , and c i are, respectively, the electric potential, Faraday's constant, relative permittivity, valence, and diffusion coefficient of ion species i in solution.J i is the total ion flux and u, , p, and  are, respectively, the velocity, density, pressure, and viscosity of water. [3,4]The specific boundary conditions in Figure 1c are as follows: applied bias voltage (line A -B), ground (line G -H), axisymmetric boundary (line A -H), bulk ion concentration (line A -B and G -H), and external pressure (line A -B and G -H).The detailed boundary conditions are summarized in Table S2 (Supporting Information).Moreover, the model being developed here takes advantage of the default particle-boundary interaction models available in COMSOL (details are summarized in Note S2, Supporting Information).

Analyte Translocation
The motion of a streptavidin-like sphere can be calculated according to, where, m a and v are, respectively, the mass and velocity of the sphere, while t time.Exerting on the sphere, F D is the drag force, F E is the electrophoretic force, and F B is the Brownian force.The force terms are expressed as follows, [3,4,27,36]  where, d a , e, Z, E, , k B , T, and Δt are the diameter of the sphere, elementary charge, charge number, electric field, generated random vector, Boltzmann constant, the absolute temperature of the fluid, and time step taken by the solver, respectively.The electroosmotic force (F EOF = 3d a u) is included as part of F D .The spheres in the simulation were released at the initial time (t = 0 s) and the path of each sphere was computed for a total length of 10 ms at an interval of 0.25 μs, the choice of Δt is detailed in Note S3 (Supporting Information) along with the physics of  in Note S4 (Supporting Information).

Electroosmotic Flow (EOF) Inside TCP
In the simulations to follow, a bias voltage was set to the upper reservoir in the range of −500 to +500 mV stepped by 100 mV, with respect to the ground set to the lower reservoir.The simulation results shown in Figure 2a represent typical EOF patterns and vortexes formed at −500 mV (left half) and +500 mV (right half).It is well-established by simulation [8] that the forcefully drifted ionic flow in the EDL inside the nanopore is the root cause for the formation of the EOF vortex (Figure 2a; Figures S2-S4, Supporting Information).For all the results, cold color (including blue) indicates the upward flow while warm color (including red) denotes the downward flow.The straight orange lines in Figure 2a mark, in the two bias conditions, where the narrowest distance from the nanopore sidewall to the boundary of the vortex is located.On the boundary lines represented by the white curves in Figure 2a, the flow velocity in the z-axis direction is zero.Inside each vortex, the EOF is in the opposite direction to that along the sidewall due to the incompressibility of electrolytes.The z-direction flow velocity along the orange lines is found in Figure 2b to first rapidly increase from zero on the nanopore sidewall, to peak at ≈1 nm from the sidewall, and then to decrease to zero on the vortex boundary.
At fixed salt concentration and surface charge density, the narrowest gap between the nanopore sidewall and the vortex, defined by the length of the orange lines, is seen in Figure 2c for a 10 nm diameter (d p = 10 nm) TCP to significantly decrease with increasing the applied bias voltage of both polarities.The results in Figure 2c conclude that the EOF vortexes move downward closer to the nanopore orifice as the bias voltage (both polarities) increases, which concurrently leads to an almost linear increase in average z-direction flow velocity of the EOF at the narrowest gap displayed in Figure 2d.These trends prevail for larger pore diameters, but with weaker EOF and more distanced vortexes from the nanopore orifice (see Figures S3 and S4, Supporting Information for the simulation data of a TCP of d p = 18 nm).

Electric Field Along the Nanopore Axis
In a nanopore-electrolyte system, the applied bias voltage can affect the transport of analytes via two distinct mechanisms, i.e., electrophoresis and electroosmosis.The detection of translocating analytes is realized by analyzing temporal drops of the ionic current through a nanopore as schematically shown in Figure 1b for such a case.The current drop is intuitively caused by the temporal interception of a translocating analyte with the most resistive section of the nanopore.At least 90% of the total resistance of the nanopore-electrolyte system results from this section confined within the so-called effective transfer length, L eff , defined as the sum of distances where the electric field falls to e −1 of its maximum on both sides of the nanopore along the central axis. [37]he highest electric field of ≈108 V m −1 is found near the orifice (z = 0) of a TCP of d p = 10 nm in Figure 2e.The electric field rapidly decreases as it moves away from the orifice resulting in L eff = ≈16 and ≈29 nm, respectively, for TCPs of d p = 10 and 18 nm (Figure 2f).The L eff of TCPs are shorter than the membrane thickness, but they are ≈2-3 times the average spatial step of the randomly moving analytes (≈9.2 nm) due to F B when adopting a timestep of 0.25 μs for calculating F B based on Equation (7).The size of F B predicted using this equation is obviously dependent on the Δt chosen, but the reported values here are found to be in very good agreement with experimental results of the root-mean-square space step of particles over time. [34]n terms of spatial steps (Table S3, Supporting Information) that are also determined by the Δt chosen, they should be appreciably smaller than the nanopore dimension and the vortex size in order to ensure sufficient details of interactions of translocating analytes with their immediate surroundings.In this regard, the 0.25 μs timestep used appears sufficient for capturing, with certainty, a trajectory and an eventual analyte translocation in a TCP.

Dominant Forces for Analyte Translocation
It is plausible that a translocation event starts from an analyte migrating from the bulk electrolyte toward the nanopore predominantly via diffusion (details in Note S5 and Figures S7 and S8 of Supporting Information).The analyte, then, in the presence of an external electric field, drifts through the nanopore thereby completing the translocation process.The breakpoint between diffusion and drift can be roughly defined by L eff in Figure 2f.The highly concentrated electric field inside L eff results in a similarly sharp force distribution peaking inside and falling by 2-3 orders of magnitude immediately outside the nanopore. [9,32]Thus, the contributing forces for the drift are F EOF and F E .Outside the nanopore including the bulk electrolyte, the analyte diffusion is determined by F B .This theoretical analysis is supported by tracking the evolutions of the forces with time and correlating them to the trajectories of five randomly selected analytes all assigned a negative charge with Z = −5.346(Note S6 and Figure S9, Supporting Information).During most of the translocation time, F B (Figure 2i) quickly alternates between ≈0.01 and 10 pN, while F EOF (Figure 2g) and F E (Figure 2h) stay below ≈1 fN until the analytes intercept the nanopore during which F EOF becomes comparable to F B and F E is appreciably greater than F B .The consequence of F E being the dominant force inside L eff to, at the final step, drive the analytes through the nanopore is a total disagreement between simulation (Note S7 and Figures S10-S13a, Supporting Information) and experiment [8] with respect to the translocation direction.
The Z value of analytes is known to be affected by the salt concentration as well as the pH of the electrolyte they are in.Employing the procedure reported in previous studies [38,39] led to the calculated effective Z of streptavidin (Protein Data Bank ID: 1SWE) to be ≈−1, −2, and −3 for 1×, 0.1×, and 0.01× PBS, respectively.The two types of streptavidin (3RY1, details in Note S6 and Figure S9, Supporting Information, and 1SWE) were calculated to carry the same Z by themselves according to the online protein calculator (https://www.protpi.ch/Calculator/ProteinTool), but the actual effective Z ought to be lower due to the Debye screening effect. [39]The screening effect has also been observed when studying protein surface charge. [40]By varying Z from −0.1, to −1, and to −5.346, the final step of analyte translocation was found to successively shift from being F EOF -dictated to F E -driven inside L eff as manifested by the distinct change of the translocation direction (Figures S12 and S13b, Supporting Information).Therefore, the consideration of the high salt concentration of 5× PBS used in the experiment [8] motivates our simplification of streptavidin as a charge-neutral analyte (Z = 0) in the remainder of the present work.

Analyte Translocation
A successful analyte translocation can be defined with the assistance of L eff (Note S5 and Figures S7 and S8, Supporting Information).Briefly, a translocation represents a complete process of an analyte initially (t = 0 s) placed in the upper (lower) reservoir traveling toward the nanopore, subsequently passing through the pore region defined by L eff , and finally landing in the lower (upper) reservoir.With this definition, the dwell time (t 2 -t 1 , Figure S8, Supporting Information) of analytes is the duration of their presence in the region confined by L eff .The mean dwell time will be discussed momentarily below.The FTE of the charge-neutral analytes (simulated three times at each bias voltage) was typically 20-50 per 10 ms for both d p = 10 and 18 nm TCPs in the range of −500 to +500 mV bias voltage (Figures S14 and S15 Supporting Information).The translocations are clearly EOF-dominated and the translocation direction is identical to that of EOF shown in Figure 2a.Those few translocations that occur in the opposite direction to the EOF can be attributed to the randomness in motion of analytes by F B and represent a negligibly small fraction compared to the number of translocated analytes by F EOF (Figures S14 and S15 of Supporting Information).They arise due to Equation (7) developed for the common "diffusive Brownian motion," but it has a limitation in describing the true trajectory of an analyte since a very short timestep (≈ 0.1 -1 μs) would be required to observe a "ballistic Brownian motion." [41,42]his incompatibility constitutes an inevitable error in order to apply F B with the other forces on the same timestep (Δt = 0.25 μs).Despite the relatively rare F B -driven translocations, a number of iterative simulation calculations were carried out to improve the statistical weight.Hence, only the translocated analytes by F EOF were registered to obtain the FTE ratio in Figure 3a,b for TCPs of d p = 10 and 18 nm, respectively.The FTE ratio for both TCPs dis-plays a similar trend converging to unity with increasing (absolute) bias voltage in the range of 100-500 mV.This phenomenon is caused by another vortex formed in the lower reservoir when a negative bias voltage is applied and increased (Figures S2 and S3, Supporting Information).At negative bias, the analyte translocation mediated by F EOF is determined to be in the direction from the bottom to the top.The increase in both dimension and flow velocity of this lower vortex with the amplitude of negative bias voltage narrows the gap between the nanopore and the vortex, thereby, making it increasingly difficult for the analytes to enter the nanopore from below.The observation of the simulated FTE ratios converging to unity at high bias voltages for both TCPs of d p = 10 and 18 nm is in qualitative agreement with the experimental data. [8]However, clear quantitative discrepancies remain to be addressed.The difficulty in fully replicating realworld experiments in simulation is an outstanding challenge.For instance, changes in the salt concentration of the electrolyte will not only alter the magnitude of F EOF but will also affect the effective net charge of the analytes and should be considered for both F EOF and F E , while the effect of partially varying viscosity within the nanopore should be carefully incorporated into the simulation as it will affect F EOF , F E , and F B .Approximating the experimental TPP using a TCP to facilitate the time-demanding simulation may have overlooked some subtle yet critical details such as the variation of viscosity with dimension (Note S1 and Figure S1, Supporting Information).The qualitative agreement is, nonetheless, considered sufficient for employing the model to perform the analysis of the forces acting on the analytes.The trajectory and corresponding forces during the dwell time of a randomly selected analyte are depicted in Figure 3c,d for the TCPs of d p = 10 and 18 nm at a fixed bias voltage of +500 mV, respectively, to elucidate the important details at the final step of a translocation (along with similar data but at −500 mV in Figure S16, Supporting Information).In both cases, F EOF plays a decisive role whereas F B also has a significant impact, and F E = 0 due to charge-neutral analytes assumed.Along the trajectory from the upper to lower boundaries of L eff , F B dominates over F EOF and F E during most of the dwell time with a random movement of the analyte inside the region between the boundaries of L eff .In other words, the analyte will largely move randomly for F B > F EOF while maintaining the directional movement by F EOF ; but it will follow the streamline of EOF flow with slight wobbling for F B < F EOF .While the magnitude of F B , |F B |, fluctuates in the range of 0.05-0.94pN (Figure 3c,d), the analyte can often reach a spot at which the magnitude of F EOF , |F EOF |, exceeds |F B |. Finally, the analyte completes the translocation by the dominant F EOF at ≈1.5 pN.This trend becomes more pronounced for larger d p (Figure 3d), as F EOF weakens rapidly upon its movement away from the sidewall (Figure 2b).
The ability to trace the trajectory of individual randomly moving analytes (Figure 4 4c,f.These statistical results support the presence of an internal vortex as the cause responsible for the increase in dwell time because the EOF inside the vortex is in the opposite direction to the EOF along the nanopore sidewall (Figure 2a).The distribution of dwell time is Gaussian for all cases but the one with d p = 10 nm at +100 mV (Figures S17-S20, Supporting Information), which leads to the extraction of mean dwell time summarized in Figure 4g,h.For both cases of d p = 10 and 18 nm, the mean dwell time tends to decrease with increasing bias voltage, conforming with the anticipated impact of F EOF whose magnitude is seen to increase with bias voltage in Figure 4a,d.At ±100 mV, the analytes experience a lower EOF impact than at higher voltages because |F EOF (max)| always lies in the range of |F B |, implying that the analytes retain the directionality by F EOF but have a more random behavior than otherwise at higher bias voltage.Hence, the analytes have a longer mean dwell time from ≈3 to ≈6 μs.Above ±100 mV, the mean dwell time decreases and it becomes ≈1 (d p = 10 nm) and ≈2 μs (d p = 18 nm) at ±500 mV.Finally, the simulated mean dwell times are also comparable to experimental results for other types of proteins, [22,23,43] implying that our model can be helpful in exploring the experimentally inaccessible details due, for instance, to bandwidth limitations as well as in developing SSNPs with new design such as TPP and BNP. [22,23,43]

Conclusion
We have developed a model based on numerical simulations to account for our earlier experimental observations of rectified protein translocations in geometrically asymmetric nanopores.
A key feature of the model is the consideration of Brownian force, in addition to electrophoretic and electroosmotic forces incorporated in numerous previous studies of ionic transport in nanopores.Extensive evaluation of various boundary conditions and modules available on the commercial simulation platform used, COMSOL Multiphysics, as well as of additional models is necessary and has led to a successful implementation of the model for numerical simulations.In the large parameter space considered, electrolyte viscosity is found to be an influential factor in analyte translocation behavior and outcome.The simulation predicts the electroosmotic force being determinant at the last moment of translocation for neutral analytes, insomuch as the electrophoretic force being determinant for charged ones.However, the pathfinding and statistical analysis of the analyte movements confirm that Brownian force on the analytes plays the dominant role all the way until the last moment of their translocation.The simulated dwell times show comparable values to the experimentally determined results for proteins with similar molecular weights.Our results are expected to be useful not only for the analysis of experimental results but also for the design of new solid-state nanopore sensors.

Figure 1 .
Figure 1.Model setup used for simulation with a truncated-conical nanopore (TCP).a) Schematic of a two-terminal sensing system based on the TCP of orifice diameter d p = 10 nm biased with an externally applied voltage set to the upper reservoir relative to the ground to the lower reservoir.b) Relative ionic current drop, ∆I/I 0 , with I 0 being the ionic current in the open pore state, caused by a spherical analyte of diameter 6 nm intercepting the TCP.Each point is taken from COMSOL simulation by successively changing the position of the analyte at +200 mV.c) Geometric and interface settings using the cylindrical coordinate (r, z) for a 2D axisymmetric system for COMSOL simulation.d) Enlarged geometric shape near the nanopore orifice from (c), where  0 is the surface charge density, h the thickness of TCP, d p is the diameter of the orifice, and  the angle of the sloped sidewall.In all simulations, the electrolyte and  are fixed to be 5x phosphate-buffered saline (PBS) and 54.7°, respectively.

Figure 2 .
Figure 2. Simulated ion-related properties for the TCP of d p = 10 nm.a) Distribution of electroosmotic flow (EOF) at ±500 mV, with the yellow lines marking the shortest distance between the vortex and the TCP sidewall.b) Distribution of EOF velocity corresponding to the yellow lines in (a) as a function of the distance from the sidewall.c) Dependence of the narrowest gap on bias voltage for d p = 10 nm.d) Variation of average EOF velocity in the narrowest gap with bias voltage for d p = 10 nm.e) Distribution of electric field for d p = 10 nm at −500 mV.f) Definition of effective transfer length, L eff , according to the distribution of electric field for d p = 10 and 18 nm.Comparison of the magnitude of g) electroosmotic force (|F EOF |), h) electrophoretic force (|F E |), and Brownian force (|F B | with timestep (∆t) of 0.25 μs) exerting on analytes with a negative charge of Z = −5.346,all at −400 mV for d p = 10 nm.Five analytes were randomly selected among the translocated ones to demonstrate the force evolution.Additional graphs labeled A, B, and C highlight the respective force immediately prior to the translocation event.

Figure 3 .
Figure 3. Frequency of translocation events (FTE) of charge-neutral (Z = 0) analytes.Simulated FTE ratio as a function of bias voltage for TCP of a) d p = 10 nm and b) d p = 18 nm.Time-dependent trajectory (upper panel) and corresponding forces (lower panel) of an analyte within its dwell time for TCP of c) d p = 10 nm and d) d p = 18 nm, both at +500 mV.

Figure 4 .
Figure 4. Statistical analysis of force and dwell time for translocated analytes (Z = 0).Comparison of the maximum force magnitude of F EOF , F EOF (max), to the magnitude ranges of F E and F B exerting on individual analytes at bias voltages from ±100 to ±500 mV for TCP of a) d p = 10 nm and d) d p = 18 nm.Distribution of dwell time for translocated analytes at bias voltages from +100 to +500 mV for b,c) d p = 10 nm and e,f) d p = 18 nm.Extracted mean dwell time as a function of bias voltage for g) d p = 10 nm and h) d p = 18 nm.The mean dwell time for the case of d p = 10 nm at +100 mV was taken from the geometric mean of the histogram, while it was determined by Gaussian fitting for the rest.
; Figures S7 and S8, Supporting Information) in combination with the definition of dwell time based on the concept of L eff allows for the exploration of translocation characteristics in depth.First, how the maximum magnitude of F EOF , |F EOF (max)|, varies with bias voltage within the dwell time is shown in Figure 4a,d for analytes translocating TCPs of d p = 10 and 18 nm, respectively.Whereas |F B | randomly fluctuates in the range of 30 fN-0.96 pN (filled area in red) regardless of the bias voltage, a substantial increase in |F EOF (max)| from ≈0.3 to ≈1.9 pN (obtained from the average value of translocated analytes) are proportionally dependent on the bias voltage.The relationship between |F B | and |F EOF (max)| is naturally associated with the distribution of dwell time, as shown in Figure 4b,c,e,f for d p = 10 and 18 nm, respectively.The larger distribution of |F EOF (max)| at positive bias than that at negative bias in Figure 4a,d leads to a wider distribution of dwell time in Figure