Precise capture and dynamic relocation of nanoparticulate biomolecules through dielectrophoretic enhancement by vertical nanogap architectures

Toward the development of surface-sensitive analytical techniques for biosensors and diagnostic biochip assays, a local integration of low-concentration target materials into the sensing region of interest is essential to improve the sensitivity and reliability of the devices. As a result, the dynamic process of sorting and accurate positioning the nanoparticulate biomolecules within pre-defined micro/nanostructures is critical, however, it remains a huge hurdle for the realization of practical surface-sensitive biosensors and biochips. A scalable, massive, and non-destructive trapping methodology based on dielectrophoretic forces is highly demanded for assembling nanoparticles and biosensing tools. Herein, we propose a vertical nanogap architecture with an electrode-insulator-electrode stack structure, facilitating the generation of strong dielectrophoretic forces at low voltages, to precisely capture and spatiotemporally manipulate nanoparticles and molecular assemblies, including lipid vesicles and amyloid-beta protofibrils/oligomers. Our vertical nanogap platform, allowing low-voltage nanoparticle captures on optical metasurface designs, provides new opportunities for constructing advanced surface-sensitive optoelectronic sensors.

As the authors mentioned in their discussion of previous studies using nanogap electrodes, this concept is not entirely new, although this study is the first one showing a vertical nanogap geometry with the flexibility to adapt a variety of surface designs over a large area, providing a platform that would enable optical detection over large areas.
I believe this work will be of great interest, since this approach show a novel fabrication design to exploit the VNEs advantages with a platform that allows for the fabrication of VNEs over large areas. Is still unclear to what degree the fabrication methods allow for consistent VNEs to be assembled over these large areas. Nonetheless, the capacity to generate such high electric field gradients with low potentials in a stable (no significant changes in temperature and solution due to direct contact with the electrodes) system, is something that a lot of researchers are looking into in order to manipulate NPs for a wide range of applications without the need of tagging or pretreating the samples of interest. The VNEs described in this work will shift the idea of conventional methods of generating electric field gradients, from horizontal electrodes, 3D carbon or coated electrodes, insulating structures, to these vertical nanogaps.
The manuscript in its present form needs to be proofread (I have included a highlighted copy of the manuscript and supplementary information), to improve its layout and of the supplementary information. There is no Introduction/Background section heather, and the figures of the supplementary Information are either miss-referenced and positioned in an order that makes the reading hard going back and forth.
The authors should include repetitions of their experiments with the different (or a subset of) NPs and provide a standard deviation or error for the normalized intensity or threshold potential. Additionally, I would be very interested in seeing the device robustness in terms of the fabrication process (i.e. evaluate the large area performance of the device by measuring FL intensity at different regions of the VNEs array).
The authors can also mention how they prepared the 200 um layer of PDMS (spin coating, fixed volume over fixed area) in the Supplementary Information.

First, the message of the paper sounds blurry. I understand that the author performed a gradual testing of their device from the largest objects to the smallest, but an extensive dedicated study focused on one of these types would have been more convincing since it would have allowed a more thorough analysis of the processes at stake, and also a more exhaustive presentation of the results. Here, the reader is drowned under a bunch of results which have no real connecting thread.
The authors agree with reviewer #1's critics that the previous manuscript was lacking of exhaustive analysis on particle dynamics as well as device performance generating ACEO despite its gradual testing of particle manipulation from the largest object to the smallest. Especially, it lacks information related to the electro-osmotic slip velocity and dielectrophoretic phenomena with respect to a variety of particulate properties (surface charge, shape etc.) and interactions with the captured particles.
As reviewer #1 pointed out, the authors found that the previous manuscript could mislead our main intention of this work; Introducing applicability of our vertically-aligned nanogap electrodes (VNE) to wide range of nanoscale biomolecules based on overall interpretation, majorly assisted by fluidic flows (ACEO) and capturing forces (DEP). Now in present version, the authors have provided additional analysis explaining how and why those two major forces interplay to relocate the particles in the level of tens of nanoscale. The authors believe and hope that overall interpretation of our VNE device's performance through numerical analysis will highlight the novelty of this work and solve the issues that reviewer #1 pointed out. Thank you again for giving us opportunity to improve quality of our manuscript and providing advices to extend depth of the work.
Please, check the additional comments addressed one by one below.

Second, the device itself (without particles) would have deserved a dedicated analysis. The rectified electro-osmotic/electro-hydrodynamic flow has to be characterised more accurately. Because the particles are not passive, I understand that micro-PIV is probably irrelevant, but an hydrodynamical numerical analysis of the cell using Comsol (Stokes mode) and Electro-osmotic slip velocity at the boundaries is indispensable before going further. I reckon that this is partially done in the supplementary materials, but this should be part of the main article. If this is not possible, it means that the format of Nature Comm is not suitable for the purpose.
The authors agree with reviewer #1's opinion that more dedicated analysis on electro-osmotic flow is necessary. The authors have conducted additional simulations to characterize ACEO-driven local slip velocities at the boundaries of two electrodes as well as resultant bulk fluidic flow across the insulator of VNE, compared to the case of horizontally-aligned nanogap electrodes (HNE). In the course of additional simulations, the authors found out that asymmetric electrode configuration of VNE generates unbalanced slip velocities at each surface of two electrodes resulting in a pair of hydrodynamic bulky flows in different magnitudes and directions. In contrast to the symmetric HNE's equivalent slip velocities and micro-vortices, our VNE provides unique electrokinetic flows with asymmetric slip velocities and resultant asymmetric micro-vortices. This offers a clue to investigate particles' hydrodynamic flow-assisted behaviour at different applying frequencies (i.e., particle mixing at 1 kHz and island-type particle relocalization at 10 kHz).  "The VNE array, which is designed to generate well-controlled flow dynamics accompanied with improved particle capture, …" TEXT ADDED: Main manuscript Page 5-6 "Together with DEP performance, ACEO-driven local slip velocities at the boundary of two electrodes and resultant bulk fluidic flow across the insulator of both VNE and HNE were explored. The simulations of hydrodynamic bulky flow over the each electrode revealed that VNE generates a pair of micro-vortices rotating in opposite directions with distinct magnitudes (Fig. 1h), in contrast to symmetric HNEs, which yield bulky flows of equivalent magnitude in an outward direction from the nanogap (Fig. 1i). As is well known, the slip velocity, uslip, arising at the electrode surface from the interaction of the E-field and EDL, drives the bulky fluidic motions. We analysed uslip to characterize ACEO behaviour as a function of both frequency and position on the VNE surface (Figs. 1j-1l). On the boundary of the VNE surface, uslip reached a maximum at the edge became decreased away from the nanogap irrespective of the applied frequency (Fig. 1j). Although the simulation results show good correspondence with those of previous studies examining ACEO-assisted slip velocity above HNEs 32,37-39 , the VNE offers unique hydrodynamic flows resulting in asymmetry in magnitudes of uslip at each micro vortex (Figs. 1k and 1l). This imbalance of a faster uslip on the top electrode than that on bottom electrode arises from the structural configuration of VNE, as previously demonstrated by the use of asymmetric horizontal electrode arrays 40-42 or non-planar ones 43-45 . Owing to this faster uslip on the top electrode, the micro-vortex on the top electrode induces larger flow streamlines than those on the bottom electrode (Fig. 1h) "In regime I, pairs of asymmetric ACEO flows circulating above each surface of VNE were simulated both in 1 kHz and 10 kHz. Because ACEO is predominant in the low-frequency range of regime I, the movements of particles are dominantly determined by ACEO rather than DEP. Within regime I, distinct particle dynamics are also expected with respect to the applied frequencies owing to the varying dominant forces; ACEO is dominant over DEP at 1 kHz (Fig. 3a), whereas they are comparable at 10 kHz (Fig. 3b), as shown through simulations (Supplementary Section 4). In contrast, the formation of the EDL layer is weakened at high frequency as the alternation of electrical signals is too rapid for ions to follow; thus, particles are majorly governed by FDEP (Figs. 3c-3e)."

The Authors newly added figures and explanation about ACEO-induced bulk fluidic micro-vortices
TEXT ADDED: Page 8 "At 1 kHz (Fig. 3f), the particles were swirling over the single unit of VNE because they are significantly influenced by ACEO flows rather than the considerably weaker DEP force. However, at 10 kHz, suspended particles are migrated to the centre of the unit cell, where flows are converged to form a point of stagnation ( Fig. 3g and Supplementary Sections 4). Unlike at 1 kHz, an application at 10 kHz generates both ACEO and DEP with comparable magnitudes. With an aid of attracting FDEP, particles are brought into proximity with the VNE electrode surface against upward ACEO flow. Then, because the ACEO-induced slip velocity (pushing toward centre) is found to be stronger than the DEP force (attracting toward edges) at the collecting surface of the VNE, effective surface streams push particles toward the stagnation point where the potential energy is the lowest (Supplementary Section 12) 51 . In contrast, at 100 kHz ( Fig. 3h) and 1 MHz (Fig. 3i) in regime II, the ring-shaped particle assembly was observed at the edge of nanogap, where the value of ∇E 2 is the largest and the potential energy is the lowest (Supplementary Section 12). Although ACEO flows still exist at 100 kHz, trapping at the stagnation deteriorates as a pair of asymmetric micro-vortices vanish into a single large flow and FDEP becomes stronger than FACEO (Supplementary Sections 4). Instead, a longrange ACEO flow constantly conveys distant particles to the VNE, and the dominant FDEP successfully snatches the moving NPs from the flow and accumulate them to the edge of VNEs at 100 kHz. "

No rigorous equation is provided for the dynamics of the particle. Technical details on the forces at stake are provided in the supplementary materials and if I get it right, the authors write the total force as the sum of a dielectro-phoretic force and an electro-osmotic, which is the drag force induced by the hydrodynamic flow, itself coming from electro-osmotic slip velocity at the walls… Hard to tell from the body of the article. Supplementary materials seem to suggest that this is the case. In brief, the author should provide a) the rigorous form of what they call the electro-osmotic force b) a comprehensive equation for the dynamics of the particle (even if the main purpose of the article is not theoretical) -is inertia relevant (probably not), what is the Stokes number for these polystyrene particles? Etc.
The authors agree with reviewer #1's concern that total force described in previous manuscript is insufficient to express dynamics of suspended particles in the water environment. It is now revised by adding Langevin equation with various external phenomena including dielectrophoresis, AC electro-osmosis, electrothermal flow, gravitation, buoyancy, Brownian motion, and interparticulate interactions. According to the calculation in Supplementary information, the authors found that the effects of electrothermal flow, gravitation, buoyancy, Brownian motion are negligible in our case. For low concentration of suspended particles floating apart each other, interparticulate interactions are also negligible except for special cases (the special cases will be discussed in 1.7 below). Thus, the exerting forces on suspended particles can be expressed in two governing phenomena; DEP and ACEO. The detailed information of each exerting forces on a single particle is expressed in the newly revised main manuscript. Besides, FACEO was defined rigorously in terms of Stokes drag force and addressed in detail together with dimensionless analysis on Stokes numbers of particles (Stk). All the details related to this issue will be addressed below. Note the accordance that the authors addressed in the order of reviewer #1's question (b) and (a) in the main manuscript.

In brief, the author should provide … b) a comprehensive equation for the dynamics of the particle (even if the main purpose of the article is not theoretical)"
In the main manuscript, the comprehensive dynamics of particles are given as between neighbouring particles (Fi,j; interaction force acting on the i-th particle owing to the j-th particles). Because Stokes drag forces are defined as F = -6πηR(up -um), where up and um are the velocities of particles and fluidic flows, respectively, FACEO and FETF exerting on a single particle can be evaluated by inserting uACEO and uETF into um, respectively 36 . Considering the simulations (FETF, FDEP, and FACEO), calculations (Fgrav, Fbuoy, and ξ(t)), and low concentration (10 ppm) condition (Fi,j) on the VNE having pattern size (L) and periodicity (P) of L = 10 μm and P = 30 μm (Fig. 1b), FDEP and FACEO dominantly determine the movement of 1-μm-diameter polystyrene (PS) particles, while the others are negligible ( Supplementary  Sections 1-5). Moreover, owing to the small Brownian timescale for relaxation of the particle (τB = m/6πηR ≈ 6 × 10 -8 s), equation (3) can be simplified into terminal up as 36 . (4) with the assumption that the initial velocity of the particle is 0."

"In brief, the author should provide a) the rigorous form of what they call the electro-osmotic force"
The rigorous form of FACEO is given and expressed in the main manuscript and supplementary information in detail. "Note that the electrohydrodynamic drag forces acting on the particle by ACEO were calculated using FACEO = -6πηR(up -uACEO), with assumptions that E-field had just been applied and that the particle was at rest prior to the AC application (up = 0) 46   For the quantitative comparison of simulation results, the electrohydrodynamic drag forces acting on the particle owing to ACEO (FACEO) and ETF (FETF) are calculated using F = -6πηR(up -um), in which up and um represents particle velocity and fluid velocity induced by ACEO (uACEO) and ETF (uETF), respectively. For calculating FACEO and FETF, up = 0 was initially assumed, indicating that electric potential has just been applied, as mentioned in the main manuscript 22 . Meanwhile, DEP-induced particle velocity (uDEP = FDEP/6πηR) is also calculated and compared with uACEO and uETF. These results (Fig. S4) show that ETF is negligible (in the order of 10 -18 N and 10 -10 m/s) compared with DEP and ACEO (in the order of 10 -18 N and 10 -10 m/s)."

"is inertia relevant (probably not), what is the Stokes number for these polystyrene particles? Etc."
The importance of Stokes numbers on dynamics of particles is discussed and calculated in Supplementary Information. Also, the results and meanings are briefly mentioned in the revised main manuscript.
TEXT ADDED: Supplementary Information Page 11 "Section 6: Stokes numbers of particles.
The Stokes number (Stk) is the dimensionless number that characterizes the particle behaviour suspended in a fluid flow, which is defined as the ratio of characteristic time of suspended particle to that of the fluid 24 Stk = 2ρpuR 2 /9ηL.
A particle with Stk >> 1 maintains its initial trajectory, without much deflection, instead of following the streamline of flow. On the contrary, for Stk << 1, the motion of the particles is dominated by the viscous force with negligible inertial effect, and therefore, is tightly coupled to the fluid streamlines. Since velocities arising from ACEO dynamics were calculated to be in the order of ~10 -4 m/s (Section 4 and Main Fig. 3), an Stk << 1 is calculated for the moving PS particles with various diameters used in our experiments. For example, Stk of PS particles with diameters of 1 μm, 200 nm, 100 nm, and 50 nm are Stk ≈ 7 × 10 -7 , 3 × 10 -8 , 7 × 10 -9 , and 2 × 10 -9 , respectively. Therefore, we conclude that suspended particles of our experiments move along the fluid streamlines, instead of detaching from the streamlines where flow abruptly changes."

TEXT MOVED/REVISED/ADDED: Page 4
Calculated Stokes numbers with physical meaning is briefly commented.
"From a dimensionless analysis of particles that exhibit small Stokes numbers (Supplementary Section 6), the viscous fluidic effect is expected to dominate over the negligible inertial effect and the motions of particles will be coupled with flow. Therefore, the movements and positions of the suspended particles under various AC conditions can be anticipated by calculating two dominant electrokinetics of DEP and ACEO."

Based on comments made in point 3, another question can be raised concerning the applicability of the so-called model to three species which are physically and chemically very different. Comments should be made on the relevancy of the model towards each type of particles. Possible local surface charges should be commented.
The authors are very pleased to improve our manuscript through adding more data which address the device performance in particle captures with different physicochemical properties. As reviewer #1 pointed out, dielectrophoretic response of particles with different physicochemical properties should be addressed for the bio-and nano-particulates in this article.
As a first step, the authors chose the appropriate dielectric models for particles with different structural characters and discussed in detail. By referring previous pioneering studies, Maxwell-Wagner-O'Konski (MWO) model (for PS particle, B. Subtillis, and Amyloid beta protofibrils/oligomers), multi-shell model (for spherical yeast cell), elongated multi-shell model (for spheroidal yeast cell), walled cell model (for vesicle), and elongated model (for Amyloid beta fibrils) were utilized based on their structural similarities. Based on these models, CM factor values were calculated as a function of applied frequency.
Taking a step further, the role of local surface charge on a sub-micron particle is also investigated together with its dielectrophoretic properties in the presence of applied voltages. Since the ionic currents and its convection within the EDL play an important role in polarization and induced dipole of nanoscale particle with its size comparable to EDL thickness, MWO theory was utilized to characterize these in terms of surface conductivity.

"Based on comments made in point 3, another question can be raised concerning the applicability of the so-called model to three species which are physically and chemically very different. Comments should be made on the relevancy of the model towards each type of particles."
Dielectric where εp * (ω) and εm * (ω) denote the complex permittivity of the particle and medium, respectively. Each complex permittivity is expressed as, in which εp, εm, σp, σm are the permittivities and conductivities of the particle and medium, respectively (j = √-1). Equations (12) and (13) indicate that the polarizing behaviour is majorly governed by the material properties of the particles and medium conditions. Additionally, the local surface charge of the particles can be evaluated by adopting the Maxwell-Wagner-O'Konski (MWO) model, which is given by 25 , where σbulk and Ksurf denote the bulk conductivity and surface conductance of the particles, respectively. Applying the parameters specified in Table S2 to equations (12)(13)(14), the real parts of the CM factor (Re[fCM(ω)]) were successfully calculated using their effective conductivities of σp = 114 μS/cm, 570 μS/cm, 1140 μS/cm, and 2280 μS/cm for particles with diameters of 1 μm, 200 nm, 100 nm, and 50 nm, respectively (Fig. S9)."  "Section 14: Trapping of yeast on VNE.

Preparation of yeast cells.
To prepare the yeast cells, Yarrowia lipolytica ATCC MYA-2613 strain was pre-inoculated into yeast synthetic complete (YSC) medium containing 6.7 g/L yeast nitrogen base (YNB), 20 g/L glucose, and complete supplement mixture (MP Biomedicals, Solon, OH, USA), and cultivated for 24 h at 28 ºC. The cells were then inoculated into a 250 ml flask containing 25 ml of the respective medium, and cultivated at 28 ºC for 144 h with orbital shaking at 200 rpm. The elongated cells were prepared by exposing spherical cells to stress condition with DI water for 3 h.
Calculation of CM factor of yeast Cells. To simplify the complexity of the biomaterial structure, the multishell model was adopted to evaluate the effective dielectric response of spherical yeast cells 29-31 . Based on the two-shell model consisting of three distinct regions of highly conductive cytoplasm, ultrathin lossy lipid membrane, and cell wall (Fig. S13a), the CM factor is calculated using homogeneous effective permittivity of the model. The CM factor is given as follows, and the parameters are given by, , However, for elongated yeast cells, the spheroidal shell model with a long axis radius of Rx and short axis radii of Ry = Rz was utilized ( Fig. S13b) 32,33 . In this model, the CM factor is given by [33][34][35] (17) in which A denotes the depolarizing factor. Assuming that the E-field is parallel to the major x-axis for the ellipsoid of Rx > Ry = Rz, A is given by 34,35 where e = [1 -(Ry/Rx) 2 ] 1/2 . While Re[fCM(ω)] in the spherical model is bounded within the range of -0.5 < Re[fCM(ω)] < 1, that of the spheroid model presents much greater Re[fCM(ω)] values, indicating that the DEP force acting on the elongated particle is much stronger than that on the spherical ones (Fig. S13c). Since the attractive DEP force is stronger in the direction of the long-axis of spheroid, most elongated yeasts cells are aligned in the radial direction (Fig. S14)."       "Calculation of CM factor of B. subtilis. The CM factor of B. subtilis was calculated using MWO theory, simplifying the complexity of the biomaterial structure into a homogeneous dielectric solid (Fig. S15). Note that surface conductance was assumed from those of virus particles 39-41 ."  "Calculation of CM factor of SUV. To evaluate induced dipole polarization of SUVs, the walled cell model of Jones, which consists of three different regions of cell interior, cell membrane and the wall, was adopted (Fig. S18a) 45 . By replacing the wall and cytoplasm of the model with a lossy lipid membrane shell of finite thickness and an aqueous medium, respectively, the complex permittivity of the dielectric model (εp * (ω)) is given by, (19) where each parameter is defined as, εim * (ω) = εimj(σim/ω), εLB * (ω) = εLBj(σLB/ω), and with their values specified in Table S5. By substituting εp * (ω) of equation (12) into that of equation (19), the CM factor of SUV was predicted as a function of the applied frequency (Fig. S18b)."       "Calculation of CM factor of Aβ42 assemblies. To achieve dielectric properties of Aβ42 fibrils, an ellipsoidal approximation was employed to simplify the high aspect ratio of the fibrous structure into an ellipsoidal model (Fig. S20a) 55 . Using equations (17) and (18), the Re[fCM(ω)] of the ellipsoidal Aβ42 fibrils (Fig. S20b) was estimated with a depolarizing factor of A = 0.0047. In this process, the surface conductance of the Aβ42 fibrils was also considered by using the effective surface conductivity of a charged ellipsoidal particle when the long-axis is parallel to the E-field, which is given by 56 , in which I1 = (1 -sinh 2 ξ0)·arcsin(1/coshξ0) + sinhξ0 and tanhξ0 = Ry/Rx.
While CM factor of the Aβ42 fibrils was calculated by ellipsoidal approximation, the spherical dielectric model together with the MWO theory was employed for Aβ42 oligomers/protofibrils (Fig. S20b). Owing to the lack of literature, we assumed the surface conductance of Aβ42 from those of biological macromolecules 57 .
Although the DEP force is much stronger in the direction parallel to the long-axes of spheroids, the captured Aβ42 fibrils are observed to be aligned along the VNE edge. This is attributed to their greater structural flexibility ( Fig. S21 and Main Figs. 5f and 5g) compared with stiff yeast cells (Fig. S14)."

TEXT REVISED: Main manuscript Page 9
Comments on dielectric model of yeast cells and B. subtilis.
"As an example of micron-scale bioparticles, the possible capturing ranges of yeast cells and Bacillus subtilis spores were theoretically evaluated using spherical-and elongated multi-shell dielectric models using MWO theory, and the trapping behaviours were experimentally confirmed under low voltage amplitudes of Vpp = 0.5 and 0.8 V, respectively (Supplementary Sections 14 and 15)."

TEXT ADDED: Main manuscript Page 9
Comments on dielectric model of PS nanoparticles.
"According to the calculation of the CM factor based on the MWO model, the values of fχ for each NP were calculated as 12, 18, and 37 MHz, respectively (Supplementary Section 10). From this analysis, a frequency of 10 kHz was selected to achieve ACEO dominance, and 100 kHz guaranteed DEP dominance."

TEXT ADDED: Main manuscript Page 10-11
Comments on dielectric model of SUVs.
"To predict the frequency dependent dielectric behaviour of SUVs, the CM factor was calculated by utilizing the walled cell model designed by Jones (Supplementary Section 17), followed by experimental corroborations."

TEXT ADDED: Main manuscript Page 11
Comments on dielectric model of Aβ42 assemblies.
"Prior to the experiments, the CM factors of Aβ42 peptides were calculated by using dielectric models depending on molecular shapes and sizes (Supplementary Section 19)."

"Possible local surface charges should be commented."
The authors explained and replied to this concern in #1.4.1 using MWO theory. The physics of local surface charges and its effect on dielectric behavior of nanoparticles are now revised and newly included in the main manuscript.

TEXT REVISED AND ADDED: Main manuscript Page 7
Detailed explanations on effects of local surface charge on dielectric properties of nanoparticles are addressed in the revised version using MWO theory.
"The direction of FDEP can also be controlled by changing the frequency because FDEP highly depends on the polarizing behaviour of particles. For example, particles move along the increasing field gradient at the condition of pDEP where Re[fCM()] > 0, whereas repulsion of particles from the region of the highest field gradient is referred to as nDEP where Re[fCM()] < 0 34 . In particular, for charged particles with sizes comparable to their EDL thickness, ionic currents and their convection within the EDL play an important role in polarization and formation of induced dipole 47 . To evaluate the dielectric response of a suspended particle in the presence of an applied E-field, the Maxwell-Wagner-O'Konski (MWO) model was adopted to consider both local surface charges and the resultant surrounding EDL. In this model of MWO dielectric dispersion, the effective conductivity of particles (σp) with a radius of R can be described as σp = σbulk + 2Ksurf/R, which is a combination of the bulk conductivity (σbulk) and surface conductance (Ksurf) of particles 48 . This implies that NPs with high surface area to volume ratio experience higher Ksurf; thus, surface conductance become dominant in the case of smaller NPs. Because MWO theory works at low conductive environments 49,50 , the effective dielectrophoretic properties of PS NPs, including CM factor, were represented using the MWO model (Supplementary Section 10)."

Again, related to the points previously listed, no trapping scenario is proposed. Most of the time, trapping results from the balance between two forces, but could also result from the localisation of the particle at a minimum of potential related to one single particular force (as it does for electrostatic trapping). A more thorough and rigorous analysis of the physics would help regarding this.
The authors strongly agree with reviewer's comment that the particle localization within the electrokinetic forces could be expressed in terms of minimum of potential energy (U). As an effort to explain localization of particle on particular region of the electrodes, the authors have conducted a calculation of the potential energy landscapes along the electrode. In this process, electrostatic potential energies under particle capturing conditions (f = 10 & 100 kHz) are calculated from inverse gradient of FDEP and FACEO, using relation of F = -∇U. The detailed description is newly represented in Supplementary Information and briefly mentioned at the revised main manuscript.

The potential energy along the surface of electrode was discussed and its values were calculated under the application of frequencies at f = 10 kHz and 100 kHz.
TEXT ADDED: Supplementary section 12 "Particles suspended in an aqueous environment move and migrate to positions at which the potential energy is minimum. As two dominant forces affect the local position of the suspended particles on the VNE surface, the total potential energy profile along the surface, U(x), can be given as U(x) = UDEP(x) + UACEO(x), where UDEP and UACEO indicate that the potential energies arise from DEP-and ACEO driving forces. Note that two dominant potential energies are achieved from the inverse gradient of FDEP and FACEO using the equations UDEP = -πεmR 3 Re[fCM(ω)]E 2 and FACEO = -∇UACEO, respectively 23,28 . Shown in Fig. S10, the minimum value of total potential energy U(x) is calculated at the centre of the VNE at f = 10 kHz (Fig. S10a), in contrast to the edges of the VNE at f = 100 kHz (Fig. S10b) "Then, because the ACEO-induced slip velocity (pushing toward centre) is found to be stronger than the DEP force (attracting toward edges) at the collecting surface of the VNE, effective surface streams push particles toward the stagnation point where the potential energy is the lowest (Supplementary Section 12) 51 .
In contrast, at 100 kHz (Fig. 3h) and 1 MHz (Fig. 3i) in regime II, the ring-shaped particle assembly was observed at the edge of nanogap, where the value of ∇E 2 is the largest and the potential energy is the lowest (Supplementary Section 12)."

Brownian motion and molecular diffusion are kept silent throughout the paper. I agree that they are probably negligible in most cases, but quantitative comments must be made on this specific point, especially if the particles' Stokes numbers are small.
The authors agree with reviewer #1 that quantitative analysis of Brownian motion and diffusion on our proposed VNE should be made for analyzing particle dynamics in aqueous environment. In this regard, dimensionless analysis on Péclet number for mass transport is carried out. Péclet number is a non-dimensional number, which is defined as a transport ratio advection-to diffusion rate. Therefore, under simultaneous interaction of external driving forces and Brownian diffusion, the study of particle transportation in a continuum can be estimated using Péclet number for mass transfer. The authors newly added comments for this issue in the main manuscript and supplementary information.

TEXT ADDED: Supplementary Information Page 10
Calculation of Péclet number (newly added). where L and u represent characteristic length of the system (in the order of ~10 μm) and fluid velocity, respectively. Since the simulated velocities were typically in the order of ~10 -4 m/s (Section 4 and Main Fig. 3), Pe of PS particles was calculated to be Pe >> 1 for most of the cases; Pe ≈ 2 × 10 3 , 4 × 10 2 , 2 × 10 2 , and 1 × 10 2 for PS particles with its diameters of 1 μm, 200 nm, 100 nm, and 50 nm, respectively. Therefore, we concluded that the effects of diffusion are negligible relative to the driving force, thereby the stochastic term ξ(t) can be neglected from the Langevin equation 23 ." "To evaluate the competitive forces and resultant dynamics of a spherical particle in an aqueous environment, the Langevin equation of particle velocity (up) was employed 35 : (3) considering DEP force (FDEP), Stokes drag force induced by ACEO flow (FACEO) and electrothermal flow (FETF), gravitational force (Fgrav), buoyant force (Fbuoy), random Brownian force (ξ(t)), and interaction between neighbouring particles (Fi,j; interaction force acting on the i-th particle owing to the j-th particles). Because Stokes drag forces are defined as F = -6πηR(up -um), where up and um are the velocities of particles and fluidic flows, respectively, FACEO and FETF exerting on a single particle can be evaluated by inserting uACEO and uETF into um, respectively 36 . Considering the simulations (FETF, FDEP, and FACEO), calculations (Fgrav, Fbuoy, and ξ(t)), and low concentration (10 ppm) condition (Fi,j) on the VNE having pattern size (L) and periodicity (P) of L = 10 μm and P = 30 μm (Fig. 1b), FDEP and FACEO dominantly determine the movement of 1-μm-diameter polystyrene (PS) particles, while the others are negligible ( Supplementary  Sections 1-5)."

Under AC fields, particles will acquire an electric dipole moment (especially at high frequency, where Debye layer cannot screen the charges), which will trigger dipole interactions and subsequent interaction forces (attractive or repulsive depending on the relative positions of each particles, but most of the time attractive). Again, no comment is made about this potentially cluster generating process.
As reviewer #1 noted, induced dipole moments of particles trigger interactive forces. Although their interactions are strong between adjacent particles, they become negligible for distant suspended particles with relatively low concentration. In this end, inter-particulate forces are considered negligible in Langevin equation, and the simulations were carried out without them. However, since inter-particulate attracting interaction occurs in-between highly populated particles, comments should be made on the experimental situations of particle trapping. The interactive clustering is commented in the main manuscript, but quantitative analysis is not conducted since it is out of the scope of the article,

Inter-particulate interactions are newly commented in Langevin equation and experimental part of the main manuscript.
TEXT ADDED: Main manuscript Page 4 "To evaluate the competitive forces and resultant dynamics of a spherical particle in an aqueous environment, the Langevin equation of particle velocity (up) was employed 35 : (3) considering DEP force (FDEP), Stokes drag force induced by ACEO flow (FACEO) and electrothermal flow (FETF), gravitational force (Fgrav), buoyant force (Fbuoy), random Brownian force (ξ(t)), and interaction between neighbouring particles (Fi,j; interaction force acting on the i-th particle owing to the j-th particles). Because Stokes drag forces are defined as F = -6πηR(up -um), where up and um are the velocities of particles and fluidic flows, respectively, FACEO and FETF exerting on a single particle can be evaluated by inserting uACEO and uETF into um, respectively 36 . Considering the simulations (FETF, FDEP, and FACEO), calculations (Fgrav, Fbuoy, and ξ(t)), and low concentration (10 ppm) condition (Fi,j) on the VNE having pattern size (L) and periodicity (P) of L = 10 μm and P = 30 μm (Fig. 1b)

There are big gaps in the literature. Articles by A. Ajdari, Castellanos and Squires, among others, should be mentioned.
TEXT MOVED/REVISED: Relationship between fluid velocity into force (Page 9 in previous version) → As a description in Supplementary Section 4 (re-ordered and revised).            TEXT ADDED: Description on CM factor Aβ42 assemblies (newly added).

The authors should include repetitions of their experiments with the different (or a subset of) NPs and provide a standard deviation or error for the normalized intensity or threshold potential. Additionally, I would be very interested in seeing the device robustness in terms of the fabrication process (i.e. evaluate the large area performance of the device by measuring FL intensity at different regions of the VNEs array).
The authors sincerely appreciate reviewer #2 on mentioning the device's reliability for evenly capturing over a large area. Since our VNE emphasizes its large-area applicability, the authors agree that the previous manuscript was insufficient to emphasize large-area performance consistency in a view point of FL uniformity. In this regard, the authors have re-analyzed the VNE performances from multiple-unit cells and provided statistical results with their intermediate process of analysis. The newly added data obtained from large-area analysis confirm a uniformity as well as massive capturing performance of the proposed VNE.

For trapping PS nanoparticles, time-lapse variation of normalized FL intensities and threshold voltages were analyzed on multiple unit cells of N = 49 and presented with their statistical results. Also for SUVs, same procedure was carried out on multiple unit cells of N = 81. The normalized FL intensity on a single unit cell was replaced into that of averaged values on multiple unit cells. The standard deviations of threshold voltages were calculated and given as the error bars in inset graphs.
Main Figures:   (Figs. 4b and 4e). As the AC voltage increased, the 100 nm particles began to be captured at Vth = 0.72 V while no 50 nm particles were monitored (Figs. 4c and 4e). Finally, the trapping of the 50 nm particles was observed for Vth = 1.24 V (Figs. 4d and 4e) while other larger particles became highly packed (see the SEM image of Supplementary Section 16). In this process, the time-lapse normalized FL intensities were analysed and averaged for a VNE unit cell number (N) of N = 49, and the standard deviations of Vth were calculated to be 0.11, 0.12, and 0.22 V for PS NPs of 200, 100, and 50 nm, respectively (Fig. 4e). Note that a glitch occurred at the moment when electrical signals changed owing to equipment limitations." TEXT REVISED/ADDED: Main manuscript Page 10.
"As the applied Vpp increased by 0.1 V every 10 s at a fixed frequency of 100 kHz, the SUVs began to be captured from Vth = 1.76 V (N = 81, standard deviation of 0.23 V), and the amount captured increased accordingly with time and voltage amplitude (Fig. 5a)."

"Additionally, I would be very interested in seeing the device robustness in terms of the fabrication process (i.e. evaluate the large area performance of the device by measuring FL intensity at different regions of the VNEs array)."
The Note that x and y are integers that represent the x-and y-coordinates of the pixel positions. Then, the distribution of the FL intensity achieved from a single cell in the m-th column and n-th row can be expressed as I(m,n)(x, y) = I(x + X(m -1), y + Y(n -1)) where 0 < x < X and 0 < y < Y ( Fig. S19a; Step 1). By carrying out image-processing and numerical analysis ( Fig. S19a; Step 2), the average FL distribution (Iavg(x, y)) was calculated from every I(m,n)(x, y), as follows .  Fig. S19a; Step 3). It is given by, .
For minimum and maximum FL intensities of 0 and 255 (arbitrary unit), Iavg values were calculated to be 66 and 57 for f = 10 and 100 kHz, respectively (Figs. S19b and S19e). After obtaining average FL intensities of each unit cell (I(m,n)) as follows ( Fig. S19a; Step 4), normalized FL intensities of each unit cell (I * (m,n)) were calculated by dividing each I(m,n) by Iavg ( Fig. S19a; Steps 5 and 6), over 20 × 20 array (Figs. S19c and S19f). The standard deviations of 0.054 (f = 10 kHz) and 0.052 (f = 100 kHz) indicate that 90 % of the unit cells exhibit variations in FL intensity within ± 10 %, demonstrating the large-area uniformity of the proposed VNE array (Figs. S19d and S19g)."

The authors can also mention how they prepared the 200 um layer of PDMS (spin coating, fixed volume over fixed area) in the Supplementary Information.
The authors agree that the process of preparing thin PDMS film is insufficient, and its fabrication method is described in more details.
The manuscript of Supplementary Information was revised to provide detailed fabrication process of PDMS film.

"The authors define a potential for a dissipative force or a set of dissipative forces (hydrodynamic interactions with the EO flow and DEP mainly). While I understand why they are doing this, they should mention that this is not a real potential (forces are not conservative). They could use the term 'pseudo-potential' for instance to stress that this potential is not defined in a rigorous way mechanically speaking."
[Author response] We appreciate the Reviewer's suggestion for concern that a set of forces suggested in the manuscript are not conservative. As reviewer #1 suggested, a terminology of 'pseudo-potential (which is also a dissipation function)' is a more rigorous and accurate representation of our interpretations to characterize the role of dissipative (non-conservative) forces acting on localization of particles. Now in the revised manuscript, the author uses 'pseudo-potential' instead of 'potential', and its relevant explanation were newly addressed both in main manuscript in briefly and supplementary information with depth explanation.
[Author actions] Explanation on introducing 'pseudo-potential' was newly included and words 'potential' were revised into 'pseudo-potential'.
TEXT ADDED/REVISED: Main manuscript Page 9 "Then, because the ACEO-induced slip velocity (pushing toward centre) is found to be stronger than the DEP force (attracting toward edges) at the collecting surface of the VNE, effective surface streams push particles toward the stagnation point where the pseudo-potential energy is the lowest (Supplementary Note 12) 51 . In contrast, at 100 kHz (Fig. 3h) and 1 MHz (Fig. 3i) in regime II, the ring-shaped particle assembly was observed at the edge of nanogap, where the value of ∇E 2 is the largest and the pseudo-potential energy is the lowest (Supplementary Section 11)." TEXT ADDED/REVISED: Supplementary information Page 18 "Supplementary Note 11: Pseudo-potential for localization of particles.
Particles suspended in an aqueous environment move and migrate to positions at which the potential energy is minimum. In this approach with non-conservative forces, the particle behaviour of pseudo-potential was introduced to interpret role of dissipative forces on localizing behaviours of particles. As two dominant forces affect the local position of the suspended particles on the VNE surface, the total pseudo-potential energy profile along the surface, U(x), can be given as U(x) = UDEP(x) + UACEO(x), where UDEP and UACEO indicate that the pseudo-potential energies arise from DEP-and ACEO driving forces. Note that two dominant pseudo-potential energies are achieved from the inverse gradient of FDEP and FACEO using the equations UDEP = -πεmR 3 Re[fCM(ω)]E 2 and FACEO = -∇UACEO, respectively 23,28 . Shown in Supplementary Fig. 10, the minimum value of U(x) is calculated at the centre of the VNE at f = 10 kHz (Supplementary Fig. 10a), in contrast to the edges of the VNE at f = 100 kHz (Supplementary Fig. 10b), respectively. This shows great consistency with the experimental demonstration shown in Main Fig. 3."

"The Langevin equation is written in the inertial full form. The authors should say early in the derivation (actually as soon as they write the equation) that they are aware that the ballistic regime is of short duration, the viscous regime being quickly reached. Perhaps I missed something and the authors already mentioned this, but I could not find it."
[Author response] We appreciate the Reviewer's comment. As reviewer #1 pointed out, transient ballistic regime of Langevin equation (short time scale of t ≪ τ) reaches into saturated viscous regime (long time scale of t ≫ τ) where τ denote relaxation time of particles. Since the τ of the particle is calculated to be much smaller than the typical time of observation, ballistic regime quickly reaches into viscous regime. For time scale of viscous regime (t ≫ τ), e -t/τ term can be negligible and the equation is saturated into right term in the equation (4). The authors newly added and revised comments for this issue in the main manuscript.
[Author actions] Comments on relaxation time and transition from ballistic regime into viscous regime were newly included in the main manuscript.
TEXT ADDED/REVISED: Main manuscript Page 5 "Since the relaxation time (τ) of the particle (τ = m/6πηR ≈ 6 × 10 -8 s) is much smaller than the typical experimental observation time (t), transient ballistic regime of short time scales (t ≪ τ) quickly reaches viscous regime of a time scale longer than τ (t ≫ τ). Thus, equation (2) with the assumption that the initial velocity of the particle is 0."

"I do not fully understand the notation ΣFij. I admit that the interactions can be reduced to pair interactions and that n-body interactions can be neglected, but the physical ingredients contained in this term are not clear. I assume that it contains electric/electrophoretic mutual interactions and also hydrodynamic interactions. Is that all? anyway, it should be specified."
[Author response] The authors thank Reviewer #1 for the comments to improve our work adding a depth of rigorous analysis. We now agree with Reviewer's opinion that notating particle interaction in a previously revised form of 'ΣFij' can cause confusions to the readers since the Langevin equation has been defined by forces acting on a single particle. In revised version, we replace the 'ΣFij' to 'Fint', which is described as a Columb interparticle interaction force to define particle-particle interaction acting on a single particle. This notation follows the reference, Lapitsky, D. S. Particle separation by alternating electric fields of quadrupole type. J. Phys.: Conf. Ser. 774, 012178 (2016).
[Author actions] The notation was revised in to 'Fint' and the following notations are also revised accordingly.
TEXT ADDED/REVISED: Main manuscript Page 4-5 "To evaluate the competitive forces and resultant dynamics of a spherical particle (with its mass of mp) in an aqueous environment, the Langevin equation of particle velocity (up) was employed 35 : considering DEP force (FDEP), Stokes drag force induced by ACEO flow (FACEO) and electrothermal flow (FETF), gravitational force (Fgrav), buoyant force (Fbuoy), interparticle force from Coulomb interaction (Fint), and random Brownian force (ξ(t)). Because Stokes drag forces are defined as F = -6πηR(up -um), where up and um are the velocities of particles and fluidic flows, respectively, FACEO and FETF exerting on a single particle can be evaluated by inserting uACEO and uETF into um, respectively 36 . Considering the simulations (FETF, FDEP, and FACEO), calculations (Fgrav, Fbuoy, and ξ(t)), and low concentration (10 ppm) condition (Fint) on the VNE having pattern size (L) and periodicity (P) of L = 10 μm and P = 30 μm (Fig. 1b), FDEP and FACEO dominantly determine the movement of 1-μm-diameter polystyrene (PS) particles, while the others are negligible (Supplementary Notes 1-5)." Main manuscript Page 9 "Note that the interparticle force component (Fint) is negligible when suspended particles are distant from each other."