Identification of DNA Bases Using Nanopores Created in Finite-Size Nanoribbons from Graphene, Phosphorene, and Silicene

The success of graphene for nanopore DNA sequencing has shown that it is possible to explore other potential single-atom and few-atom thick layers of elemental 2D materials beyond graphene (e.g., phosphorene and silicene). Using density functional theory, we studied the interaction of DNA bases with nanopores created in finite-size nanoribbons from graphene, phosphorene, and silicene. We observe that binding energies of DNA bases using nanopores from phosphorene and silicene are generally smaller compared to graphene. The band gaps of phosphorene and silicene are significantly altered due to interaction with DNA bases compared to graphene. Our findings show that phosphorene and silicene are promising alternatives to graphene for DNA base detection using advanced detection principles such as transverse tunneling current measurement.


I. INTRODUCTION
Research into methods of deoxyribonucleic acid (DNA) sequencing has been improving with the goal of providing fast and cheap sequencing of longer strands as accurately as possible. One possible way of improving the process is utilizing solid-state devices made from two-dimensional (2D) materials. Graphene remains the most widely studied 2D material for DNA sequencing. Graphene nanopores and nanogaps have been successfully used for DNA sequencing. 1,2 While the single-layer nature of graphene provides the optimal thickness (0.34 nm) for single-base resolution, 3 the major hindrance is the hydrophobic nature of graphene's surface. Because of the strong π-π interactions between graphene and the DNA, bases stick to its surface 4,5 leading to a significant reduction in translocation speed due to pore clogging. 6 Furthermore, the coexistence of different bases on the surface and pore makes single-base discrimination difficult. 6 Another issue is the problem of orientational fluctuations of nucleobases during DNA translocation through a graphene nanopore. This can give rise to overlapping current contributions from different bases. 7 It has been shown that the ionic blockade signal shows noise for DNA translocation through a single-layer graphene nanopore. 5 The origin of this noise has been attributed to the atomic thickness of the pore. It is notable that a nanopore in a three-layer graphite structure, which has a thickness ~1 nm, shows a better signal-to-noise ratio (SNR). 5,8 The lack of a band gap in pristine graphene makes it undesirable for use in electronic-based detection modalities such as tunneling current or field-effect transistors (FETs). 9 Despite the challenges of graphene, graphene's success for nanopore DNA sequencing has shown that it is possible to explore other potential single-and few-atom thick layers of 2D materials beyond graphene, and also that these materials can exhibit fascinating and technologically useful properties for DNA base detection that are superior to those of graphene. One such material is molybdenum disulfide (MoS2). Among the large family of 2D transition metal dichalcogenides (TMDs), MoS2 is the most widely studied for sequencing applications, mainly due to the ease of fabrication of MoS2 devices. [10][11][12] Several theoretical and experimental studies have demonstrated sequencing using single-layer nanopores and nanoribbons from MoS2. 9,10,[13][14][15] These studies reveal that MoS2 performs better than graphene. For instance, improved SNR, non-stickiness of DNA to MoS2 surface, and the presence of an intrinsic band gap makes it suitable for use in advanced sequencing devices such as FETs or the nanochannel device. 9,10 Another 2D material that has been explored for DNA sequencing is hexagonal boron nitride (hBN). hBN is less hydrophobic than graphene which can minimize the hydrophobic interaction that impedes the DNA translocation through the constructed nanopores. Moreover, the thickness of hBN is comparable to the spacing between nucleotides (0.32 -0.52 nm) in single-stranded DNA (ssDNA). 16 It also shows other advantages over graphene in terms of its insulating property in high-ionic-strength solution and fewer defects made during the manufacturing process. 17 Experimental and theoretical studies show that hBN can reach spatial resolutions for DNA sequencing close to those of graphene. [18][19] The fundamental properties of DNA nucleobases and hBN sheets remain unchanged upon adsorption, which suggests its promising application for DNA research. 20 The above consideration shows that hBN is a nice alternative to graphene and could be used in combination with graphene for advanced functional devices for DNA sequencing with higher SNR. 19 While tremendous success has been reported with the use of 2D materials (e.g., graphene, hBN, MoS2, and WS2) for DNA sequence, there are still fundamental limitations related to the high translocation speed during DNA sequencing that has to be resolved. 9,10,21 At this stage, it is not clear what solid-state nanopore material will be able to meet the challenges for single-base resolution. Thus, it is crucial to continue to carry out explorative studies to identify new nanopore materials that could potentially emerge as the best candidate materials. Two such candidate materials are phosphorene 3,22 and silicene. 11,23 The puckered honeycomb structure of phosphorene and buckled honeycomb structure of silicene makes these materials thicker than graphene, hence superior materials for single-base DNA sequencing in terms of signal-to-noise ratio intensity while maintaining monolayer properties. 9 Furthermore, phosphorene and silicene share many of the remarkable properties of graphene including high carrier mobility and tunable optical properties. [23][24][25][26] Both phosphorene and silicene have direct band gaps that are dependent on thickness, making them ideal for advanced base detection principles that probe modulations in electronic and optoelectronic properties. 23,27 Moreover, the hydrophilicity and biocompatibility of phosphorene makes it ideal for applications in biosensing. 28,29 DNA sequencing using nanoribbons from phosphorene was recently demonstrated using periodic density functional theory (DFT) calculations. 29,30 These studies were based on probing modulations in electronic properties due to translocation of DNA bases through a nanogap in phosphorene 29 or modulations caused by physisorption of bases onto the surface of phosphorene nanoribbons. 30 In this article, using DFT calculations, we studied the interaction of DNA bases with nanopores created in finite-size nanoribbons from graphene, phosphorene, and silicene. Due to the small size of DNA bases (~0.6 -0.8 nm), the bases interact with only a small section of the nanoribbon, hence using a finite-size model is appropriate for capturing the interaction between bases and 2D membrane materials. Furthermore, by using a finite-size model, our system is modeled as a molecular system, which does not require a periodic DFT calculation. We observe that binding energies of DNA bases using nanopores from phosphorene and silicene are similar, and generally smaller compared to graphene. This shows that minimal sticking of DNA bases to pore is expected for phosphorene and silicene devices compared to graphene. 4,5 Furthermore, nanopores from phosphorene and silicene show a characteristic change in the density of states for each base. The band gap of phosphorene and silicene are significantly altered due to interaction with DNA bases compared to graphene. Our findings show that phosphorene and silicene are promising alternatives to graphene for DNA base detection using advanced detection principles such as transverse tunneling current measurement. 31

FIG. 1.
Schematic of nanopore sensing device. Translocation of DNA base through pore causes variations in in-plane current. The in-plane current flows from left to right, that is perpendicular to the DNA backbone during DNA traversal.
Our model is based on detecting variations in the in-plane current 7,31 through a nanopore created in a finite-size nanoribbon membrane due to traversal of a DNA base as shown in Figure 1. First, we performed DFT studies for the individual DNA bases to determine the maximum width of each base. This was determined to be 0.754 nm for guanine (G), 0.648 nm for adenine (A), 0.560 nm for cytosine (C), and 0.590 nm for thymine (T). Then nanopores (1 -1.2 nm) were created in each of the three membrane materials, that is graphene, phosphorene, and silicene. Our pore sizes are in the 1 -2 nm range used in previous studies. 7,31 Such a pore size is good enough to allow nucleotides (size ∼ 1.0 to 1.2 nm) to pass through. 29 However, for the detection of double-stranded DNA (dsDNA), a pore diameter of 2.3 nm has been shown to be the smallest and most efficient pore size. 4,9 Nucleotides serve as monomeric units of DNA and consist of a base, a sugar, and a phosphate group. Since the sugar and phosphate groups are the same for each nucleotide, we shall consider only the interaction of the four target bases (G, A, C, and T) with nanopores from graphene, silicene, and phosphorene in this study. This approximation has been used in previous computational studies. 7,9 We will refer to our finite-size armchair nanoribbon systems using the following abbreviations: GNP (graphene nanopore), PNP (phosphorene nanopore) and SNP (silicene nanopore). Due to the small size of DNA bases, the bases interact with only a small section of the nanoribbon, hence using a finite-size model is appropriate for capturing the interaction of bases and 2D membrane materials. Furthermore, by using a finite-size model, our system is approximated as a molecular system, which does not require a periodic DFT calculation.
For each membrane material, nanopores were created by deleting central atoms. Hydrogen passivation 9,30 was used to terminate the dangling bonds at the edges and pore. Structural minimization was performed for the individual bases and membranes separately, and then for the combined system of nanopore and bases. All structural relaxation calculations were performed at the B3LYP level of theory using the 6-31G (d, p) basis set, with a force convergence cutoff of 0.02 eV/Å. 30 The combination of the B3LYP hybrid-GGA (generalized gradient approximation) functional with the 6-31G (d, p) basis set provides a reasonable basis set size for performing geometry optimization calculations. All calculations were performed using the GAUSSIAN 16 software package. 32 Computational resources were provided by the University of Central Oklahoma Buddy Supercomputing Center. 33 To evaluate the potential of phosphorene and silicene for DNA base detection, we computed the electronic structure changes induced by the presence of DNA bases inside the nanopore using DFT. Three evaluation metrics, namely, energy band gap, binding energy, and electronic density of states (DOS) were considered. As graphene is the most widely studied 2D material for DNA sequencing, our results for phosphorene and silicene were benchmarked against graphene. The energy band gap is given by where is the energy of the Highest Occupied Molecular Orbital (HOMO) state, and the energy of the Lowest Unoccupied Molecular Orbital (LUMO) state. The binding energy for the different systems was calculated as follows Here, is the binding energy for the system, + is the total energy of combined system (nanopore and base), is the energy of the nanopore, and the energy of a single DNA base. In a real device, the tunneling current signal is computed numerically from the integrated DOS. 7,34 This suggest that a modulation in the DOS of the pristine nanopore system should produce a modulation in the tunneling current. In this study, instead of performing the full calculation of the tunneling current which requires the non-equilibrium Green's Function formalism, 7,29,30,34 we shall focus on studying variations in the DOS of the pristine nanopore due to interaction with DNA bases. 9 The DOS was computed using the Lorentzian function as 35 where E is the energy, are the molecular energy levels obtained from the DFT calculation, and Γ is the linewidth of the Lorentzian spectral function. In our calculations, we used Γ = 25 meV as this corresponds to the linewidth at room temperature. 35

III.1 Graphene Nanopore (GNP).
Our finite-size GNP system has 15 dimer lines along the zigzag direction and 4 translational periods along the armchair direction. The optimized length along the zigzag direction is 1.91 nm, and the optimized length along the armchair direction is 1.80 nm. The pore diameter is approximately 1.07 nm. Table I shows the band gap and binding energies for each system. The band gap of pristine GNP is 0.221 eV. Table I shows that the band gap of GNP is not significantly affected by the presence of DNA bases. Each base slightly opens the band gap of the pristine GNP. Base T most significantly affects GNP and gives rise to an increase in band gap of 0.01 eV. The change in band gap is highest in T and lowest in G and follows the order T > C = A > G. The optimized geometry for each GNP system is shown in Figure 2 along with the corresponding binding energies.
Binding energies for our GNP calculation follow the order C >A > G ~ T. The variation of binding energy and energy gap for each base is shown in  The small variations in energy band gap and DOS of GNP due to interaction with DNA bases means that GNP is not suitable for detecting DNA bases using the in-plane transverse current method. 9,31 Also the binding energies of bases with GNP is relatively large. This indicates that the individual bases can easily stick to the surface of the graphene nanopore. Sticking to the surface can cause significant issues with base detection such as reduction in translocation speed due to pore clogging. 6

FIG. 4.
Density of States of pristine GNP and pristine GNP + base systems. EF is the Fermi energy of pristine GNP.

III.2 Phosphorene Nanopore (PNP).
Our finite-size PNP system has 16 dimer lines along the zigzag direction and 4 translational periods along the armchair direction. The optimized length along the zigzag direction is 2.48 nm, and the optimized length along the armchair direction is 1.20 nm. The pore diameter is approximately 1.03 nm. Table II shows the band gap and binding energies for each system.
Unlike the GNP, the band gap of each PNP system is significantly affected by the presence of bases. The pristine band gap for the PNP is 3.070 eV. We remark here that the band gap for the finite-size PNP system is much higher than the ~0.87 eV reported for an infinitely long phosphorene nanoribbon of identical width 30 . This difference is due to the finite size of our system, and the fact that the band gap of 0.87 eV was determined for a phosphorene nanoribbon with no pore.
From Table II, base A most significantly affects PNP and gives rise to a reduction in band gap by 0.281 eV compared to pristine PNP. Unlike other bases, base G slightly opens the band gap of PNP by 0.013 eV. The change in band gap is highest in A and lowest in G and follows the order A > T > C > G. The optimized geometry for each PNP system is shown in Figure 5 along with the corresponding binding energies for the different bases. The binding energy is a key parameter that describes the strength of the interaction between the DNA base and PNP. In GNP, the binding energies follow the sequence C > A > G > T. For PNP, the binding energies follow the sequence C > G > A > T. The variation of binding energy and energy gap for each base is shown in Figure 6.  The large variations in band gap energy and DOS, and small binding energies of PNP due to interaction with DNA bases makes PNP a better material compared to GNP for DNA base detection using the inplane transverse current method. 31

III.3 Silicene Nanopore (SNP).
Our finite-size SNP system has 13 dimer lines along the zigzag direction and 4 translational periods along the armchair direction. The optimized length along the zigzag direction is 2.57 nm, and the optimized length along the armchair direction is 2.58 nm. The pore has a diameter of approximately 1.20 nm. Due to the buckled structure of silicene, creating a perfectly circular pore was very challenging. Our calculated Si-Si bond length after a geometry optimization calculation was in the range 0.228 nm to 0.236 nm, which agrees nicely with the experimental value of 0.228 nm. 36,37 Furthermore, our calculated out of plane height for the buckled structure was found to be 0.045 nm, which also agrees nicely with the experimental value of 0.044 nm. 38

FIG. 7.
Density of States of pristine PNP and pristine PNP + base systems. EF is the Fermi energy of pristine PNP. Table III shows the band gap and binding energies for each system. The band gap of the pristine SNP system was calculated to be 2.219 eV. Base A most significantly affects SNP and gives rise to a reduction in band gap by 0.151 eV compared to pristine SNP. All bases give rise to a reduction in the band gap for SNP. The change in band gap is highest in A and lowest in T and follows the order A > C > G > T. The optimized geometry for SNP and DNA bases is shown in Figure 8. In the GNP, the binding energy for each base showed C > A > G > T. For PNP, the binding energy followed the order C > G  From Figure 10, we observe that the change in band gap and DOS is unique for each base. Also, the change in band gap is large for the SNP device than the GNR device, again confirming the hypothesis that silicene might be a superior sensing material for transverse current detection compared to graphene.

III.4 Comparing GNP, PNP, and SNP.
This research sought to explore the properties of finite-size nanoribbons from graphene, phosphorene, and silicene as DNA base detectors using the nanopore method. The pore size for GNP (1.07 nm) is comparable to the pore size for PNP (1.03 nm). The pore size for SNP was slightly larger (1.20 nm). Despite the slight differences, we can compare the general trends of the modulation in electronic properties due to interaction with different DNA bases. Table IV shows the binding energy and energy gap for the different materials. We observe from Table IV that binding energies of SNP and PNP are comparable, and lower than binding energies in GNP. To further quantify the changes in band gap, we calculated the change in band gap as the difference in band gap of the pristine system and base, and the pristine system without DNA base. Table V shows the change in band gap for each material. It shows that the modulations in energy gap for PNP and SNP are significant compared to GNP. Finally, Figure 11 shows the binding energy comparison for each of the different membrane materials in the presence of each base. The lower binding energies of DNA bases in PNP and SNP suggest minimal sticking of DNA bases to pore for these materials. Hence PNP and SNP are expected to resolve the issue of pore stickiness, which is a major challenge for DNA translocation through graphene nanopore. 4,9,31 11. Visualization of differences in binding energy (eV) for GNP, PNP, and SNP.

IV. DISCUSSION
In this study, ions and solvating water molecules were not included in our idealized simulation. 7,9,29,30 Other studies have shown that solution do play an important role. We would expect solution to change the magnitude of the binding energies and band gap, but the trend will not be changed. [39][40][41] For the orientation of bases, we focused only on the orientation with bases parallel to plane of nanoribbons (θ = 0). Orientation does affect the binding energy and band gap of the various systems. For example, the binding energy and band gap of GNP for different orientations of DNA bases are shown in Table VI and Table VII, respectively. We observe from Tables VI and VII that the modulation in binding energy and band gap of GNP for different orientations of DNA bases is not very significant. Similar results are expected for PNP and SNP. Therefore, one would expect the trends in Table V and Figure 11 to be the same for different orientations of the DNA bases.   Our finite-size model is a simplification of real systems. As a result, our values of binding energies obtained using Gaussian local basis set is expected to be different from bulk (infinitely long nanoribbons) values that are usually obtained with plane wave basis set. Binding energies calculated using Eq. (2) may be subjected to basis set superposition error (BSSE). For the relatively decent size of basis set and choice of model chemistry used in our calculations, the BSSE is expected to be small, and hence can only affect the magnitude of binding energies but not the general trend. To further support this argument, we report our recently obtained results for adsorption calculations for DNA bases on finite-size nanoribbons using the same basis set and level of theory that was used for the nanopore calculations (see Section II), and then compared our finite-size calculations with results that were recently published in the literature for infinitely long nanoribbons. 30 For simplicity, we present results for only the finite-size phosphorene nanoribbon (PNR). The finite-size PNR has 9 dimer lines along the zigzag direction and 3 translational periods along the armchair direction. Figure 12 shows the optimized geometry for each base, along with corresponding binding energies (absolute values), and equilibrium adsorption heights (all bases were placed 3.0 Å above surface of nanoribbon before geometry optimization). The binding energies follow the order G > A > C > T. This same order was observed for phosphorene nanoribbons using periodic plane wave calculations. 30 The binding energies (0.169 to 0.330 eV) for the finite-size model are lower than those obtained using periodic calculations (1.75 to 2.32 eV), but the order G > A > C > T is the same for both finite-size and periodic systems. 30 The discrepancy in binding energy magnitudes is due to the difference in supercell size used for the periodic calculation. Our calculated band gap of pristine PNR is 3.038 eV. This value is higher than the value obtained using periodic calculations (~ 1.0 eV) for same width PNR. 30 This discrepancy is due to the finite-size approximation used in our model. The variation of binding energy and energy gap for each base is shown in Figure 13. Each base opens the band gap of the pristine PNR. Base G most significantly affects PNR and gives rise to a decrease in band gap of 0.358 eV. The change in band gap is highest in G and lowest in C and follows the order G > A > C ~ T. Interestingly, this order is also observed for interaction of DNA bases with infinitely long nanoribbons from phosphorene, 30 graphene, 42 and MoS2. 9 Thus, we conclude that our choice of basis set and model chemistry appropriately captures the qualitative aspects of the systems investigated. The inclusion of BSSE correction could therefore lead to a quantitative change in binding energy for finite-size systems, but the trend is expected to remain unchanged.

V. CONCLUSION
In summary, using DFT, we studied the interaction of DNA basis with nanopores created in finite-size nanoribbons from graphene, phosphorene, and silicene. We observe that binding energies of DNA bases using nanopores from phosphorene (0.2 -0.4 eV) and silicene (0.2 to 0.3 eV) are similar, and generally smaller compared to graphene (0.8 -1.1 eV). This shows that minimal sticking of DNA bases to pore is expected for phosphorene and silicene devices. The binding energy sequence was determined as C > A > G > T for graphene, C > G > A > T for phosphorene, and C > G > T > A for silicene. Furthermore, nanopores from phosphorene and silicene show a characteristic change in the density of states for each base. The band gap of phosphorene (ΔEgap ~ 13 − 281 meV) and silicene (ΔEgap ~ 43 − 151 meV) are significantly altered due to interaction with DNA bases compared to graphene (ΔEgap ~ 7 − 10 meV). Our findings show that phosphorene and silicene are promising alternatives to graphene for DNA base detection using advanced detection principles such as transverse tunneling current measurement.