Robust signatures in the current-voltage characteristics of DNA molecules oriented between two graphene nanoribbon electrodes

In this work we numerically calculate the electric current through three kinds of DNA sequences (telomeric, \lambda-DNA, and p53-DNA) described by different heuristic models. A bias voltage is applied between two zig-zag edged graphene contacts attached to the DNA segments, while a gate terminal modulates the conductance of the molecule. The calculation of current is performed by integrating the transmission function (calculated using the lattice Green's function) over the range of energies allowed by the chemical potentials. We show that a telomeric DNA sequence, when treated as a quantum wire in the fully coherent low-temperature regime, works as an excellent semiconductor. Clear steps are apparent in the current-voltage curves of telomeric sequences and are present independent of lengths and sequence initialisation at the contacts. The current-voltage curves suggest the existence of stepped structures independent of length and sequencing initialisation at the contacts. We also find that the molecule-electrode coupling can drastically influence the magnitude of the current. The difference between telomeric DNA and other DNA, such as \lambda-DNA and DNA for the tumour suppressor p53, is particularly visible in the length dependence of the current.

In this work we numerically calculate the electric current through three kinds of DNA sequences (telomeric, λ-DNA, and p53-DNA) described by different heuristic models. A bias voltage is applied between two zig-zag edged graphene contacts attached to the DNA segments, while a gate terminal modulates the conductance of the molecule. The calculation of current is performed by integrating the transmission function (calculated using the lattice Green's function) over the range of energies allowed by the chemical potentials. We show that a telomeric DNA sequence, when treated as a quantum wire in the fully coherent low-temperature regime, works as an excellent semiconductor. Clear steps are apparent in the current-voltage curves of telomeric sequences and are present independent of lengths and sequence initialisation at the contacts. independent of length and sequencing initialisation at the contacts. We also find that the molecule-electrode coupling can drastically influence the magnitude of the current. The difference between telomeric DNA and other DNA, such as λ-DNA and DNA for the tumour suppressor p53, is particularly visible in the length dependence of the current. Following the publication of the seminal work of Fink and Schönenberger on the electrical conduction of DNA strands 1 and, shortly afterwards, a single molecule version by Porath, Bezryadin, Vries and Dekker 2 , the possibilities of electronic transport in DNA have stimulated a large body of work for more than a decade. 3 This is of course driven by the exciting fundamental characteristics of DNA as a biological structure with self-recognition and self-assembly properties. However, even after all this work, we are still far removed from a proper understanding of charge transport (CT) in DNA. Experimental work continues to remain very challenging and publicationswhile reporting tantalising hints of surprisingly large currents -have only recently begun to probe CT properties in DNA beyond the magnitude of the supported current. [4][5][6] Similarly, the multitude of theoretical studies has yet to agree on which conduction mechanism and model to choose for a proper coarse-grained description of CT in DNA -or to find the means of attacking the complexities of larger DNA strands via powerful ab initio methods.
Some of the problem stems from the fact that DNA is not just a single molecule, but rather describes the whole set of DNA strands made possible by stringing the nucleotide bases Adenine (A), Cytosine (C), Guanine (G) and Thymine (T) into the classical double helix structure. Consequently, it is difficult to compare the results of publications when the DNA sequences used are different, or, as sometimes happens, not fully specified. Furthermore, the used DNA sequences range from simple pe-riodic poly(dG)-poly(dC) variants -and of course their AT counterparts -to viral, bacterial and eukaryotic DNA and beyond into artificially generated sequences.
Here, we intend to contribute to the discussion by putting forward a naturally occurring DNA sequence as a particularly well suited test bed for DNA experiments as well as theoretical studies: telomeric DNA. In mammals, it is a Guanine rich sequence in which the pattern TTAGGG is repeated over thousands of bases. Its length is known to vary widely between species and individuals and it essentially has a buffer function at the beginning and end of DNA strands for eukaryote cells. 7 It therefore combines the advantages of a simple, periodic structure with the richness and biological function of a real DNA sequence. In addition, we shall study not a single, but rather 5 different coarse-grained, tight-binding models of DNA. In this way, and in the absence of a preferred model as discussed above, we shall not be over concerned with quantitative differences between the models, but rather aim to elucidate their qualitative agreements.
Some of the most exciting solutions to the wellestablished contact problem in nano-devices is based on the use of carbon allotropes such as carbon nanotubes (CNT). 4,5 These CNT then allow for good coupling to macroscopic gold contacts. Even simpler should be the use of graphene flakes or nano-ribbons as potential contact material. In fact, the basic building block of a DNA device might very well be a single graphene sheet into which gaps are fabricated by lithographic techniques. It is these devices which we will take as our starting point here (cp. Fig. 1). 1. (color online) Schematic representation of a laddertype model for electronic transport along DNA between two semi-infinite graphene nanoribbon source (S) and drain (D) contacts as indicated. A third (gate) terminal modulates the conductance of the molecule. The nucleotide base pairs sequence is indicated by small (pyrimidines, red, yellow) and large (purines, blue and grey) circles representing the four possible effective nucleotides, T, C, A, and G, respectively. The sugar-phosphate backbone is given by green circles and possible electronic pathways are shown as lines.

II. THEORETICAL FORMULATION
A schematic figure of a molecule of DNA coupled to graphene electrodes is shown in Fig. 1. The source and drain electrodes are identified with their respective Fermi levels µ S and µ D , the effect is taken into account using self-energy functions. 8 A gate terminal modulates the conductance of the molecule. We focus on coherent transport and hence compute G r , the retarded Green's function between the electrodes S and D via the recursive lattice Green's function technique. 9 This then allows us to calculate the density of state, transmission spectrum and the current voltage characteristics.

A. The model Hamiltonians
The study of nanoelectronic applications based on macromolecules as building blocks is a compromise between the model complexity and the length scale that the chosen model can handle. Therefore, the use of heuristic tight-binding models has grown in interest in the last decade. 10 In spite of the simplified electronic structure, very long systems can be studied with these models in contradistinction to ab initio approaches. Hence, we study here the transport properties of DNA according to 5 such heuristic models. For sake of clarity, the 5 models are briefly described as follows. First we address a single strand chain with 4 different sites, representing the nucleotides (1L). 11 This simplest possible model is followed by a single strand of bases attached to sites emulating the sugar-phosphate backbones (1L+BB or "fish-bone" model). 12 Next in complexity appears a double strand chain of nucleotides (2L), 13 and a model where the said double strand of base pairs is attached to backbones (2L+BB) 14 and, finally, a double strand of base pairs with backbones and additional diagonal inter chain hoppings among base pairs (2L+BB+D). 15 Instead of repeating the mathematical definitions of the corresponding Hamiltonians as given in the above cited literature, let us concentrate on the model 2L+BB+D as the most complicated one. The other 4 models can then be reconstructed from it by suitable choice of parameters. In Fig. 1, we show a schematic representation of 2L+BB+D. The model contains parameters for the description of the backbone as well as longitudinal, transverse and diagonal base to base hopping terms which take into account the size difference between pyrimidines (C,T) and purines (A,G). The Hamiltonian can be written as l |l, α(σ) l, σ|) + h.c., (1) where t → l,j is the hopping at base pair l along the strand starting from 5' (j = 1) and 3' (j = 2) ends, t 1ց2 l and t 2ր1 l denote the diagonal hopping and t ↓ l the hopping perpendicular from 5' down to 3' at l. The sum over σ in (1) marks the connection to the sugar-phosphate backbone as in Fig. 1. Last, the hermitian conjugate indicates the hopping terms associated with t ← l,j = t → l−1,j , t 1ւ2 l . In addition, α(σ) = 1 (2) for σ =↑ (↓) and ǫ lα and ǫ The different terms are better appreciated by referring to Fig. 1. The onsite energies ǫ lα are taken to be the effective primary ionization energies of the base nucleotides, i.e. ε A = 8.24 eV, ε T = 9.14 eV, ε C = 8.87 eV and ε G = 7.75 eV. In this work, we consider the backbone energy to be given as average of the energies of the base nucleotides, i.e. ε 1 eV for purine-to-purine, 0.01 eV for purine-pyrimidine and 0.001 eV for pyrimidine-to-pyrimidine. Perpendicular couplings to the backbone sites are t [σ] l = 0.7 eV, and the perpendicular hopping across the hydrogen bond in a base pair is reduced to t ↓ l = 0.005 eV. From previous discussions leading to these choices of parameters as well as the influence of the environment on the charge migration properties of the models, we refer the reader to the existing literature. 14,17-19 Although we use rather simple Hamiltonians to describe the molecule, we believe that the qualitative physics of CT in the molecule is captured. This is because both the molecular energy levels and the wave functions closely resemble those calculated from the much more sophisticated ab initio theory. 17 Nevertheless, we emphasise that the choice of the tight binding parameters is far from uniquely determined, being a rather controversial issue, and several parameter sets have been proposed in the literature. 20

B. Green's function techniques
The transmission probability T (E) between the electrodes can be evaluated by where G r is the retarded Green's function of the system which can be found from 8 Here 1 1 is the identity matrix, U is the gate potential and Σ S , Σ D are the self-energies for source (left) and drain (right) contacts (cp. Fig. 1), respectively. These self-energies are calculated as usual from g, the Green's function of the electrode, 8 and the coupling τ between the DNA molecule and electrode (cp. Fig. 1), i.e. S,D = τ † S,D g S,D τ S,D . The Green's function g is calculated using a recursive technique. 21 The electrode-molecule coupling τ is determined by the geometry of the chemical bond. 22 We use a constant coupling τ = 0.35 eV of similar magnitude as the inter-chain DNA couplings. 2 We emphasise that our results remain robust for small changes in this parameter. The anti-Hermitian part of the self-energy is known as the broadening function and related to the lifetimes Γ S,D = i Σ S,D − Σ † S,D of an electron in a molecular eigenstate. The Hückel approach 23 predicts that the Fermi energy E F is closer to the highest occupied molecular orbital (HOMO) level. In this work the Fermi energy of the Dirac cones in undoped graphene is aligned to the mid gap (∼ 8.1 eV) of the DNA model shown in Fig. 1. Consequently, the Fermi levels in source and drain are µ S = E F and µ D = E F + eV D , respectively. The gate potential U , which allows to model the charging effects in the presence of bias, 24 can be expressed as U = eV g assuming that there is no large charge redistribution by applying bias between the electrodes.
Given H, S,D , the chemical potentials µ S,D in the electrodes and the gate potential U , we obtain the density of states DOS from G r as and the current is given as usual via In the low- where Θ is the step function. So, the electronic low-T current can be expressed as 8 The results shown in Fig. 2 are the calculated I-V characteristics of three different types of DNA sequences for short DNA strands of length L = 30 bps. The sequences considered are (i) the eukaryotic telomer based on repeats of the TTAGGG sequence as discussed in in section I, 7,25 (ii) a random subsequence of bacteriophage λ-DNA 26,27 and (iii) a random subsequence of the DNA strand of the p53 gene. 28,29 The λ-DNA sequence has been studied previously as a typical example of a biological DNA sequence. It contains 48502 base pairs with 25% of A, 24% of C, 28% of G, and 23% of T strung together in a non-periodic sequence. Similarly, p53-DNA has 20303 base pairs and is part of an important regulatory mechanism in humans. [30][31][32] B. Model dependence of the current-voltage characteristics The main results for the I-V characteristics depicted in Fig. 2 can be summarized in a few general trends. First, the inclusion of backbones opens a gap between the HOMO and the lowest unoccupied molecular orbital (LUMO), clearly revealed in the corresponding I-V curves showing threshold voltages. 33 Here it is worth recalling that we assume a Fermi energy close to half the energy gap, by tuning a gate voltage of V g = 9.8 eV for all models (see appendix A). For sufficiently high applied voltages the entire band will be scanned across the Fermi energy, leading to current saturation, as can be observed for almost all models for the voltage range depicted in Fig. 2 and further illustrated in the appendix. Next, a non-periodic sequence of base pairs, such as in λ-and p53-DNA, can drastically suppress the currents, even for such a short system of L = 30 bps, as one sees comparing Fig. 2, for disordered λ-DNA and p53-DNA, respectively, to the telomer sequence. This strong suppression is model dependent: single strands do not show a large supression at such short lengths, while double strand models already exhibit a six orders of magnitude drop in current. Furthermore, the double stranded chain 2L+BB+D, i.e. with diagonal hopping and backbones, shows higher currents than the other double stranded models, partially recovering single strand values. Last, and important in what follows, there are steps appearing in some I-V curves, which are due to resonances in the transmission probabilities. Such resonances are more relevant and robust in the telomeric sequences, due to the split of each band in several sub-bands. 14 Resonance effects are still present in some cases for disordered short sequences, but they do not last for longer systems, as can be systematically appreciated in Fig. 3. Hence we see that the DNA model proposed in (1) seems to interpo- late between the extreme cases of the simple 1L models, retains the gap due to the backbone, but has the smaller current values, that are more in line with the experimental results. In what follows, all the results are for the model (1), namely 2L+BB+D, although it is the most complicated one of those considered previously.

C. Length dependence of the current-voltage characteristics
One of the most striking differences between the telomeric and the non-telomeric sequences is revealed by comparing short and longer sequences. Fig. 3 shows the I-V characteristic for all three DNA strands, considering in each case four different strand lengths. The degree of non-periodicity in the sequences severely affects the I-V characteristics. First, there is the expected drop in the current due to non-periodicity (acting like an effective disorder), and, the robustness of the step-like structures in the I-V curves for the telomeric strands at any L.
Most prominent, however, is the rapid drop in the current when increasing L for the non-telomeric sequences. For 300 bps (104 nm), λ-DNA and p53-DNA, which are nonperiodic, support practically no current, even in the pA scale, while the telomeric sequence retains a maximum saturation current of around 0.6 µA. Fig. 4 shows the current of the three sequences as a function of length at fixed saturation current V D = 1.5 eV. For very short lengths up to L ≈ 10 bps, the three cases can not be distinguished. However, for longer strands, the telomeric sequence shows a saturation current of the order of several hundreds of nA almost constant in L, whereas for λ-DNA and p53-DNA the current decreases exponentially as the electronic states are strongly localized. It should be noted here that this result is robust and independent from the heuristic 2L+BB+D model, i.e. it happens in all the models considered. It should therefore be advantageous to consider the naturally occurring telomeres of DNA sequences as prime candidates when looking for good conductivity in a DNA strand -in nanotechnology applications but also in DNA-related in-vivo processes where charge migration might be important. 18

IV. TRANSPORT CHARACTERISTICS OF TELOMERIC DNA
In order to establish the features of the electronic transport which are intrinsic to the telomeric sequence, the role of the possible variations in coupling to the graphene contacts has to be studied. Three main aspects are considered: (i) the dependence on which base pair of the telomeric sequence is connected to the graphene, (ii) whether contacts are only made to a single strand of the DNA or to both, 3' and 5', ends, and (iii) what happens when more than a single telomer is contacted in parallel. Let us recall that for the results presented up to now, we always assumed that the transport starts at the T of the TTAGGG sequence, that the coupling τ holds for both 3' and 5' ends and that only a single telomer has been contacted.

A. Starting the telomer
In Fig. 5, we show I-V curves for telomeric sequences starting at different base pair positions, i.e. in addition to the periodic repeat of TTAGGG, we also have TAGGGT, AGGGTT, . . . , GTTAGG. In this way, we hope to capture the possibilities which might arise when the contacts are made simply by putting telomers on top of graphene sheets. 34 The DOS of the TTAGGG system, as well as T (E), is also shown in Fig. 5. The results are for L = 300 bps in all cases. The DOS shows a peak at E F related to the nanoribbon localized state at this energy. The overall envelope of the DOS is shaped by the graphene band structure and the band splitting of the DNA system. The superposed sharp peaks are related to the resonances in T (E), which show a huge variation in intensity. Therefore, some peaks in the DOS are related to transmission peaks that fall below the scale of the figure. The irregular step-like structure in the I-V characteristics reflects the fragmented structure of T (E). Notice that the apparent gap in the I-V characteristics does not necessarily correspond to the gap in the DOS. The reason is that many of the states close to the band gap have very low transmission probabilities (related to highly localized resonances at different parts of the rather long sequences), therefore they do not contribute to transport.
When inspecting the I-V characteristics for the different starting positions, one notices that each one shows steps at the same energy, although the steps show variable amplitudes. We therefore find that, while changing the start of the sequence at the contact can alter the height of the transmission probability for a given energy, it does not change the position of the steps in the I-V curves. As we had argued in section III, it is exactly the existence of these pronounced steps which makes the telomeric sequence stand out from other DNA sequences. Our result of the section now show that this finding remains robust with respect to what actual base pair is chosen to contact to the graphene sheets.

B. Contacting the telomer at 3', 5' or both
In the literature, 4,5 it is known that the choice of contacting single 3' and 5' ends only or both 3' and 5' ends can alter the current response. This is also observed in our 2L+BB+D model -the class of 1L models can of course not capture this behaviour -as we show in Fig.  6. In the figure, we have contacted the graphene electrodes via (i) only the direct 3'-to-5' or 5'-to-3' ends at source and drain, (ii) only the cross 3'-to-3' or 5'-to-5' ends and (iii) both 3' and 5' ends. Besides the trivial effect of dropping the current by reducing the number of channels, as readily observed by comparing the black curve to the others, one also sees that all cases shows a significant modification of the sub band transmission intensities in the energy range of interest around the Fermi energy. We interpret the difference between 5' and 3' results in Fig. 6 as being due to the sequence of nitrogenous bases. In the source 5' terminal, the onsite energy for the thymines (9.14 eV) is larger than the onsite energy for the adenine (8.24 eV) in the 3' terminal. So, a contact via the A in 3' is preferable. In addition, direct transmission as in 3'-to-5' can make use of the larger hopping strengths compared to the diagonal or perpendicular couplings. We note, however, that this does not fully explain the sequence of current strengths, e.g. it is not immediately clear why the source 5' to drain 5' current is larger than the source 5' to drain 3' current. Nevertheless, the present results strongly suggest that the step structure of the I-V characteristics, although intrinsic to the telomeric sequence DOS, are filtered and selected by the choice of how the strands are coupled to the contacts. Note that we still observe that the (voltage) position of the steps remains fairly robust. C. More than a single telomer in parallel The advantage of the DNA device outlined in section I and shown in Fig. 1 is the macroscopic size of the graphene contacts together with a small, nano-sized gap to be manufactured into the graphene sheet. It is then clear that it will be hard to guarantee that only a single DNA strand lies across the gap. Hence we expect that one might occasionally encounter a situation where more than one DNA strand is being contacted. As we shall show here, this seems indeed beneficial, i.e. we find that the situation of telomeric DNA strands contacted in parallel enhances the size of the current.
In Fig. 7, we show the I-V characteristics for a single sequence, together with the curves for 5 and 10 sequences in parallel, all of them L = 301 bps long, and separated each by 40 C-C distances on the graphene contacts. We start the first sequence with TTAGGG. . . and hence end after 301 bps as . . . TTAGGGT. Then we use the remainder of the telomer as the start sequence of the next strand, i.e. TAGGGT. . . AGGGTT. The end of the 5th strand is then . . . GTTAGG and the 10th strand ends as . . . GGGTTA. In this way, we have different starting and end parts of the telomer at the source and drain contacts, respectively, for parallel sequences. From Fig. 5, we can see that the step structure is indeed enhanced and not washed out by adding many strands, even if their assembly in an actual experiment is not base pair exact. It should also be noticed that the currents are not a simple linear addition of the N = 1 case, since interference ef- fects appear, due to the different sequences coupled to the same contacts. Also, the position of the steps between plateaus remains again very robust.

V. CONCLUSIONS
The complexity of the DNA systems does not allow yet a definitive conclusion to be drawn on the mechanism(s) leading to charge transport in telomeric or other sequences. Nevertheless, the consistency among the results found here for different heuristic models, widely discussed in the literature, suggests that important qualitative physical and chemical aspects have been captured in the present approach. Telomeric sequences contacted to graphene sheets may show high currents as well as robust plateau structures in the I-V characteristics. This plateau structure is enhanced and not washed out by adding parallel telomeric chains, irrespective of the starting points defined at the interface of each sequence with the contact. Most importantly, the drastic difference in the lengths dependence of the current for telomeric DNA as compared to other DNA strands such as the λ-and p53-DNA shown in Fig. 4 will allow for better comparison and interpretation of experimental data. We emphasise that the inclusion of some external disorder, be it in the fidelity of the sequence or the environmental conditions, does not drastically alter these differences. 14 An important current rectification behaviour can also be seen due to the differences in the transmission probabilities of occupied and unoccupied levels. A device including a gate electrode could further tune the rectification possibilities addressed here.
The devices investigated in the present workbased on the naturally occurring, physiochemically stable telomers -should provide ideal benchmark situations for systematic investigations of DNA electronic transport, as well as for development of DNA based molecular nanoelectronic applications using graphene as contacts. In particular, the I-V characteristics presented here show promise for a wide range of interesting nanocircuitry on complex patterns of hollowed graphene sheets bridged by telomeric DNA sequences. ACKNOWLEDGMENTS C.J.P. acknowledges financial support from FAPESP (Brazil) and the hospitality of the Centre for Scientific Computing at the University of Warwick (U.K) were most of the work has been developed. P.A.S. acknowledges partial support from CNPq (Brazil).
Appendix A: Varying the gate voltage Some features in the I-V characteristics depend on the gate voltage, V g , applied to the system and are worth briefly discussing here. The gate voltage will shift the DOS modifying drastically the transport characteristics. Keeping in mind that the main results shown in the present paper are for V g = 9.8 eV (bottom panel of Fig.  8 inset, showing the transmission probability bands superposed to the density of states) a situation where the equilibrium Fermi energy lies closer to the HOMO (top of the valence band in semiconductor terminology), one observes that negative V D show higher currents than positive V g in Fig. 8. This asymmetry is due to the lower transmission probability for the occupied molecular or- The contacts to the DNA devices are semi-infinite graphene nanoribbons as depicted in Fig. 1. Graphene is a single layer of carbon atoms packed in a honeycomb lattice, constituted by two sub lattices. Finite pieces of graphene are limited by two different kinds of edges (or, more generally, a mixing of both kinds): armchair-like edges, where the atoms at the edges alternately belong to both sub lattices and a zig-zag edge, in which all the edge atoms belong to the same sub-lattice. The graphene contacts considered here are semi-infinite nanoribbons along the zig-zag direction, hence the DNA strands are connected to an armchair edge (perpendicular to the zig-zag direction) of width W . The thinnest nanoribbons considered are W = 50 wide, in units of C-C bond length of 0.142nm. Bearing in mind that DNA chains exist in many different conformations that an average measure between 2.2nm to 2.6nm wide, the connection of an increasing number of telomeric DNA double strands lead to the use of contacts, up to W = 500, in units of C-C bond lengths, in the case of 10 telomeric DNA strands in parallel. Fig. 9 shows the energy bands of grapheme nanorribons along the zig-zag direction for two different widths, W = 25 and W = 50. As can be seen here, from the electronic properties point of view, zig-zag like nanoribbons are metallic with the valence and conduction band touching at two points as expected (armchair edged nanoribbons may be semiconductor or metallic, depending on the width). 35,36 The flat bands at the Fermi energy are localized states at the semi-infinite zig-zag edges, giving rise to the peak in the DOS at E = E F in the inset of Fig. 5, but showing no influence on the transport properties of the contacted DNA strands. Note also that we have obtained qualitatively similar results -in terms of differences between telomers and λ-and p53-DNA as well as the lengths dependence of current magnitudes and the robustness of current steps -when using simple square lattice contacts.