*Department of Molecular and Cellular Biology, Harvard University, Cambridge
Department of Chemistry, University of California, Davis
![]() |
Abstract |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Introduction |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
In a recent study, the rate constants and possible tunneling pathways for intramolecular ET between two [4Fe-4S] clusters in a series of native and genetically engineered Clostridia ferredoxins were reported (Kyritis et al. 1997
). This ferredoxin is a small (55 amino acids) monomeric protein, which contains two redox centers[4Fe-4S] cubane clusterseach ligated by four cysteine sulfurs. The lengths of the cubane edges are roughly 2.3 Å; the shortest distance between two atoms on different cubanes is 8.5 Å, whereas the shortest distance between two sulfur atoms ligating different cubanes (sulfur on Cys14 and that on Cys43) is 6.2 Å. It was found that the structural changes introduced around and between the two cubane clusters via site-directed mutagenesis had no major affect on the ET rate constant. In the same study, it was also found that the ET rate constant was independent of temperature in the range of 283310 K.
In this paper, we determine the electron tunneling pathways and calculate rate constants in the recently recrystallized Clostridium acidurici ferredoxin (Dauter et al. 1997
[The structure used for this paper was obtain from the Protein DataBank (www.rcsb.org): PDB code 2FDN]) as well as a mutant, using the rigorous quantum mechanical method of tunneling currents. We find that there is a main tunneling pathway running directly from Cys14 to Cys43 and a weaker pathway running from Cys14 via Ile23 to Cys18, whereas other amino acids do not play a significant role in the electron tunneling. In silico mutation of Ile23 to valine did not significantly affect the ET rate constant, in agreement with the experimental findings. This is because of the prominent Cys14 to Cys43 pathway.
We have also performed a molecular evolutionary study of this protein; it was found that Ile23 is a highly variable amino acid compared with the donor- or acceptor-ligating cysteines. Sequence analysis was performed on this ferredoxin, homologous sequences were found in species representing all the three domains (Woese 1982
), including metazoans and members of the deepest part of the tree of life (Lazcano and Miller 1996
). As this protein predates the prokaryote and eukaryote divergence date of around 2 billion years (Nealson and Conrad 1999
) and is ubiquitous because of its role in important biological functions, we propose that the evolutionary constraints have been placed on the cysteine amino acids not only for the purpose of setting the redox potentials of the cubanes for which these amino acids serve as ligands but also as electron tunneling mediators. These results not only rationalize experiment but also add support to the hypothesis that there exists evolutionarily conserved electron tunneling pathways in biological ET reactions (Medvedev, Daizadeh, and Stuchebrukhov 2000)
.
![]() |
Materials and Methods |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
|
Here we describe the electronic structure of the system within the extended Hückel approximation (Yates 1978
; Whangbo et al. 1987
). Extended Hückel parameters are listed in table 1
; the diagonal elements of the Hamiltonian of the system are the empirically determined ionization potentials. The nondiagonal matrix elements are given by the Wolfsberg-Helmholtz equation:
|
|
|
The transfer matrix element is calculated with the method of tunneling currents (Stuchebrukhov 1996, 1996a, 1997
). In this method, the matrix of interatomic tunneling currents is defined by
|
|
The geometric structure used for these calculations was the original crystallographic structure as determined from the experiment. In certain cases, such as conformational or energetic studies, it is advisable to equilibrate such structures by molecular dynamics or molecular mechanics simulations when the protein coordinates experience appreciable differences from those of the crystal structure. Thus, taking into account the theoretical arguments and the fact that in these systems the closest distance between the donor and acceptor complex is within 9 Å with rigid substructures including covalently bound cysteine mediated donor and acceptor complexes, no equilibration was performed. Lastly, we should mention that, as the results will show, the calculated matrix element agrees well, within the approximations used in our calculations, with that of the experimental values, thus adding an additional level of confidence in our approach.
Lastly, we note that quantum theory and simulations of the ET rate constant have shown that the matrix element becomes a sensitive function of the atomic coordinates between the donor and acceptor when the distance between the donor and acceptor complexes becomes appreciably long (Daizadeh, Medvedev, and Stuchebrukhov 1997a
). For the average intramolecular ET protein, with donor and acceptor distance around 14 Å (Page et al. 1999
), such coordinate corrections to the calculation of the matrix element are not relevant. Indeed, the discovery that the matrix element is a function of coordinate configurations was a recent theoretical finding from our group (Daizadeh, Medvedev, and Stuchebrukhov 1997a
; Daizadeh, Gehlen, and Stuchebrukhov 1997b
). Because of the close proximity between the donor and acceptor cubanes and the closeness of our calculations with that of the experimental values for the rate constant, these additional quantum mechanical simulations were not performed.
In order to calculate the donor and acceptor states, we construct the donor (acceptor) complex by taking all eight atoms of the cubane and S, Cß, H1ß, H2ßatoms on the four ligating cysteine amino acids: Cys8 (Cys18), Cys11 (Cys37), Cys14 (Cys40), and Cys47 (Cys43). The C
atoms on these amino acids are treated as hydrogen atoms and included in the donor and acceptor complexes, in order to eliminate the influence of dangling bonds on the zeroth-order donor and acceptor states. The coordinates of these hydrogen atoms are the same as those of C
atoms in the crystallographic structure.
Molecular Evolutionary Study
The entire C. acidurici primary sequence was used in the search against the nonredundant, highly curated SwissProt database using the web-version of BLAST (Madden, Tatusov, and Zhang 1996
[The web service can be found at http://www.ncbi.nlm.nih.gov/BLAST. The SwissProt database provided by NCBI was posted on December 2, 1999 and contained 82,258 sequences [Bairoch A., and R. Apweiler. 1997. Nucl. Acids Res. 25:3136.] The default BLAST parameters were used in the calculations: gap-penalties: Existence: 11, Extension: 1, with a Lambda ratio of 0.85; a low complexity filter and the BLOSUM62 scoring matrix were also used.]); 80 pairwise alignments with E-value scores less than 0.001 were harvested from the output file (the original sequence was also retrieved). The set of pairwise alignments was then multiply aligned with MULTALIGN (Corpet 1988
[Version 5.3.2 was used with the following parameters: scoring matrix = Blosum62, gap penalty = -1, gap length = -1, exremity gap penalty = no penalty, scoring method = absolute, conservation for uppercase letter in consensus = 90 for lowercase letter in consensus = 50 via the world wide web at: http://npsa-pbil.ibcp.fr/cgi-bin/npsa_automat.pl?page=npsa_multalin.html]). Because of the small number of residues that comprise this protein (55 amino acids), the sequence was shuffled, changing the order of the amino acid randomly without disturbing the overall amino acid frequency, and used to search the database. No sequences below an E-value of 0.001 were retrieved from these randomized sequences, thus adding confidence that the sequences retrieved in the search procedure using the C. acidurici primary sequence were indeed biological homologues. (In fact, in many searches no sequences were found with E-values less than 10.0; in several, E-values more than 2.0 were retrieved.)
To compare the amino acids implicated in the ET pathways by the theory summarized above and by experiment to those by sequence analysis, we developed the following heuristic approach: From the multiple alignment, the number of amino acid differences with respect to the amino acid sequence obtained from the C. acidurici PDB file (Dauter et al. 1997
) were determined columnwise. These corresponding numbers of mutations were then normalized with respect to the number of sequences in the multiple alignments (80) (see caption of fig. 4
for complete details). These numbers were then used as a basis for a coloring scheme for an image of the three-dimensional coordinates of the molecule (MSI User Guide 1995
[Figures 1 and 2
were produced by using the charge color submodule in the InsightII program with the following criterion: redness of atom p = - 2.5 x
q|Jp·q|/|TDA where intensity of color is maximum when redness is less than or equal to 0.5. Figure 5
was produced in a similar way where blueness of an amino acid is proportional to the degree of the conservation shown in fig. 4
]) to qualitatively compare the results from the theory.
|
|
|
![]() |
Results and Discussion |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
However, the SOMO of the reduced complex and the next two highest orbitals were not considered donor or acceptor states in our calculations. Electron density in these states is localized mainly on the cubane iron (40%45% of Mulliken population) and sulfur atoms (50%55% of Mulliken population), with a small contribution from the cysteine ligand atoms. However, density functional theory (DFT) calculations performed on this system (Mouesca et al. 1994
) show that both Mulliken and electrostatic potential charges of cysteine sulfurs change appreciably (by 10%15%) on each sulfur atom upon reduction of the oxidized cluster, whereas charges of iron atoms change less significantly. This shows that the donor and acceptor orbitals in our calculations should have a significant electron density on the cysteine sulfur atoms. The orbitals which satisfy this criterion are HOMO and HOMO -1 of the oxidized complex. They have 45%55% of the electron density on the cysteine sulfurs, mainly on their p orbitals, 30%35% of electron density on the iron atoms, mainly on their d orbitals, and a small contribution from other atoms. This is not in exact agreement with the results of the DFT calculations, but the discrepancy of populations of the cubane atoms is expected to affect the value of the matrix element less significantly than that of the populations of the cysteine sulfurs. The reason is that if the population on the cysteine ligands is adequately described, then the degree of delocalization of the donor and acceptor wave functions, which is crucial for the value of matrix element, is correct.
Matrix Element and Tunneling Currents
The picture of tunneling currents in the case when both donor and acceptor states are represented by the HOMO of the oxidized complex is presented in figures 1 and 2
. Tunneling currents in the case when donor or acceptor state, or both, is represented by the HOMO -1 of the oxidized complex are similar. The larger intensity of red in the figures corresponds to a larger amount of total current (the sum of incoming and outgoing currents) passing through the atom. In both the figures, the atoms with total current larger than the whole current in the system are shown with maximal color intensity. All protein atoms and directions of major interatomic currents are shown in figure 1
. In figure 2 only those amino acids with observable currents are shown; here the interatomic currents are presented in greater detail. Interatomic currents that are larger than 30% of the whole current in this system are shown with solid black arrows, whereas some of the weaker interatomic currents are depicted with dotted black arrows. Currents inside donor and acceptor complexes, except for those between atoms of the cysteine ligands, are not shown.
|
There is a strong and complicated system of circular currents involving Cys14, Cys43, and to a smaller extent the methyl group of Ile23. Another strong system of circular currents involves Cys18, Ile23, and some of the atoms of the acceptor cubane. On the other hand, Ala44 and Pro52, together with Cys47 (a donor cubane ligand) are involved in a weak system of circular currents and hence are not important for the ET. Current from the sulfur of Cys47, which has a significant amount of electron density on it flows mainly to the donor cubane and then to the acceptor via the main pathways. There are also strong circular currents inside donor and acceptor cubanes which are not shown.
The calculated value of the matrix element is between 3.7 and 6.7 cm-1. The driving force of this reaction is small, -6 and 15 meV between the two [4Fe-4S] cubane groups in C. acidurici and in Clostridium pasteuranium, respectively (Kyritis et al. 1997
). Calculated values of the matrix element correspond to the range of the rate of ET (assuming
G0 = 0) from 2.2 x 105 s-1 to 7.4 x 105 s-1 with
= 1 eV and from 4 x 107 s-1 to 1.3 x 108 s-1 with
= 0.5 eV (at room temperature). This range of values of the matrix element is in reasonable agreement with the experimentally defined rate of transfer of 6.5 x 106 s-1 (Kyritis et al. 1997
), given the crude nature of the underlying extended Hückel parameterization used in the model.
The experimental mutant, I23V, in the homologous C. pasteurianum ferredoxin did not affect the value of the ET rate significantly, there was only a marginal decrease to 5.3 x 106 s-1 from 6 x 106 s-1 (Kyritis et al. 1997
). We performed an in silico site-directed mutagenic experiment, computationally transforming Ile23 to a valine using the program SCRWL (Bower, Cohen, and Dunbrack 1997
[The SCRWL program can be obtained from: http://www.fccc.edu/research/labs/dunbrack/scwrl/]) and calculated the matrix element for the C. acidurici system. The SCRWL program was used because in our previous calculations on cytochrome c oxidase (Medvedev, Daizadeh, and Stuchebrukhov 2000)
and photolyase (Cheung et al. 1999
), the program, which replaces side chains onto the protein backbone based on predefined parameters defined from a backbone-dependent rotamer library, is a computationally nonintensive alternative to the computational procedure of replacement of amino acids followed by minimization of conformational energies, while producing the same qualitative results with respect to our electronic structure calculations. Upon this substitution, the matrix element did not change by more than a factor of 1.8; this again is in good agreement with the experiment. Examination of the tunneling currents shows that in this case almost all the current flows from Cys14 to Cys43. This shows that Ile23 is not important for the ET, but may be important for protein stability (Kyritis et al. 1997
).
We have also performed calculations using the SOMO of the reduced cluster and two next orbitals as donor and acceptor states. As expected, the matrix element was too small (not more than 0.3 cm-1) to account for the observed rate of transfer, assuming that the reorganization energy is within a usual range, which for proteins is typically taken to be between 0.5 and 1 eV. However, the amino acids participating in the electron tunneling were found to be the samecysteine ligands of cubanes, Ile23, Ala44, and Pro52. This observation adds further confidence to our results.
Sequence Analysis
Here, we present the results of the sequence analysis study of the C. acidurici ferredoxin. The homologous sequences of this ferredoxin belong to species from the three domains, which include the eukaryotes (10 sequences), bacteria (54 sequences), and the archaea (16 sequences) and span the full gamut of evolutionary distances from mammals to the green sulfur bacteria. The 10 eukaryotic sequences were composed of three Metazoa (two mammals, one nematode), three Viridiplantae, one Fungi, one Rhodophyta, one Reclinomonas, and one Entomoebidae. The 54 bacterial sequences were composed of 32 Proteobacteria, 18 Firmicutes, 3 Green-sulfur bacteria, and 1 Thermo/Deinococcus group. The 16 archaeal sequences constituted 15 Euryarchaeota and 1 Crenarchaeota. Thus the protein is an ancient protein predating the prokaryotic and eukaryotic divergence around 2 billion years ago.
The proteins, which contain regions of homology with the 55 amino acid sequence from C. acidurici ferredoxin, are associated with ET, but from a wide range of biological functions in natural systems; from an examination of the annotations for the sequences extracted from the database, we find the following distribution using the annotation support provided within SwissProt documentation: 34 ferredoxin (24 ferredoxin, 5 FDI, 4 FDII, 1 seven-iron ferredoxin fragment), 8 NADH-ubiquinone oxidoreductase 9 hydrogenases/hydrogenlyases (3 coenzyme F420 hydrogenase gamma subunits [8-hydroxy-5-deazaflavin-reducing]), 2 periplasmic hydrogenase large subunit, 1 formate hydrogenlyase subunit 6 [hydrogenase-3 component F], 2 hydrogenase-4 [including component A and H], [1 NADH-dehydrogenase I chain I], 2 DMSO reductase iron-sulfur subunit anaerobic dimethyl sulfoxide reductase chain B), and 1 RNA polymerase (DNA-directed RNA polymerase subunit D). The remaining 25 sequences correspond to sequences annotated with keywords, such as probable, putative, ferredoxin-like, or hypothetical. A complete list of accession numbers as well as the details for the multiple alignments is shown in figure 3 .
|
Because of the ubiquity and age of this protein, as well as the presence of sufficient experimental data, this protein system makes an ideal candidate to study a speculation concerning the role of evolution in the formation of tunneling pathways. A comparison of the amino acids implicated in the tunneling process by the method of tunneling currents with that of the multiple alignment data reveals an interesting correlation. In figure 4 , we show the (normalized) number of mutations of amino acids which are different from those found in C. acidurici ferredoxin, and in figure 5 , we have represented this information graphically on an image of the X-ray crystallographic structure of the ferredoxin molecule. We find that the cysteine amino acids are the most conserved part of the protein (shown in white), whereas Ile23, even though its position in the three-dimensional structure is directly between donor and acceptor cubanes, is highly variable (shown in blue). Notice that in addition to the cysteines, prolines in positions 19 and 48 are also highly conserved; the cysteines and prolines are exposed to solution and may serve as regions for protein-protein interactions. (The isopolymorphic alanincs in positions 22 and 51 are located within 4 Å of the cubane cysteine ligands and may serve to stabilize the protein).
Evolutionary substitutions of Ilc at position 23 include the hydrophobic amino acids: 9 leucines, 5 phenylalanines, 3 valines, 2 methionines, and 1 histidine. Lastly and interestingly, roughly six homologous sequences did not have a cysteine at position 11 (Cys14, Cys18, Cys37, Cys43, and Cys17 showed no variability; the substitutions at position 8 and 40 were of the order of one or two and may be because of errors in the database). It is known that homologous proteins with a deficient cysteine ligand correspond to a substitution of [4Fe-4S] with the [3Fe-4S] cofactor (Betrand et al. 1995
). Even though this substitution at position 11 will change the reduction potential, we argue that as Cys11 and Cys13 are present and it is highly likely that they are roughly the same distance apart (because of the high level of homology as it is known that sequence homologies on the order of 20%30% yield almost identical three-dimensional structures [Holmes 1999]), the ET tunneling pathway will be the same.
Although our results do not provide a logical proof that conservation of Cys14 and Cys43 is exclusively because of their involvement in the tunneling process (It can be argued that these amino acids are conserved along with the other five cysteines as ligands for regulating the redox potentials. Recall, however, that one of the ligands, Cys11 is less conserved evolutionarily than the other respective amino acids [see fig. 4
]), the conservation of amino acids making up tunneling paths is interesting and is in line with the hypothesis of evolutionary conservation of tunneling paths. Further analyses including domain-level sweeping of the database in search for other related ET is currently under investigation. This paper along with an analysis in cytochrome c oxidase (Medvedev, Daizadeh, and Stuchebrukhov 2000)
demonstrates that the union of rigorous quantum chemical theories based on crystallographic information, experiment and molecular evolution may serve together to clarify the nature of the biological ET tunneling process.
![]() |
Note added in proof |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Acknowledgements |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() |
Footnotes |
---|
Present address: Nautilus Biotech, Evry, France
Keywords: electron transfer
evolution
Clostridia
ferredoxin
sequence analysis
Address for correspondence and reprints: Alexei A. Stuchebrukhov, Department of Chemistry, University of California, Davis, California 95616. stuchebr{at}chem.ucdavis.edu
![]() |
References |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Betrand P., O. Mbarki, M. Asso, L. Blanchard, F. Guerlesquin, M. Tegoni, 1995 Biochemistry 34:1107011079.
Bower M. J., F. E. Cohen, R. I. Dunbrack, 1997 Prediction of protein side-chain rotamers from a backbone-dependent rotamer librarya new homology modeling tool J. Mol. Biol 267:1268-1282[ISI][Medline]
Cheung M. S., I. Daizadeh, A. A. Stuchebrukhov, P. F. Heelis, 1999 Pathways of electron transfer in Escherichia coli DNA photolyase: Trp306 to FADH Biophys. J 76:1241-1249
Corpet F., 1988 Multalin multiple sequence alignment with hierarchical clustering Nucleic Acids Res 16:10881-10890[Abstract]
Daizadeh I., J. N. Gehlen, A. A. Stuchebrukhov, 1997b. Calculation of electronic tunneling matrix element in proteins: comparison of exact and approximate one-electron methods for Ru-modified azurin J. Chem. Phys 106:5658-5666[ISI]
Daizadeh I., E. S. Medvedev, A. A. Stuchebrukhov, 1997a. Effect of protein dynamics on biological electron transfer Proc. Natl. Acad. Sci. USA 94:2708-3703.
Dauter Z., K. S. Wilson, L. C. Sicker, J. Meyer, J.-M. Moulis, 1997 Structure of Clostridium Acidurici ferredoxin. Detailed geometry of [4Fe-4S] clusters in a Protein Biochemistry 36:16065.[ISI][Medline]
Duéc E. D., E. Fanchon, J. Vicat, L. C. Sicker, J. Meyer, J.-M. Moulis, 1994 Refined crystal structure of the 2[4FE-4S] ferredoxin from Clostridium acidurici at 1.84 angstrom resolution J. Mol. Biol 243:683-695[ISI][Medline]
Homes K. C., 1999 Structural biology Philos. Trans. R. Soc. Lond. B 354:1977-1984[ISI][Medline]
Huber J. G., J. M. Moulis, J. Gaillard, 1996 Use of H-1 longitudinal relaxation times in the solution structure of paramagnetic proteins, application to [4Fe-4S] proteins Biochemistry 35:12705-12711[ISI][Medline]
Krauss N., W. D. Schubert, O. Klukas, P. Fromme, H. T. Witt, W. Saenger, 1996 Photosystem I at 4-angstrom resolution represents the first structural model of a joint photosynthetic reaction centre and core antenna, system Nat. Struct. Biol 3:965-973[ISI][Medline]
Kuznetsov A. M., 1995 Charge transfer in physics, chemistry, and biology: physical mechanisms of elementary processes and an introduction to the theory Gordon and Breach Publishers, New York
Kyritis P., J. G. Huber, I. Quinkal, J. Gaillard, J.-M. Moulis, 1997 Intramolecular electron transfer between [4Fe-4S] clusters studied by proton magnetic resonance spectroscopy Biochemistry 36:7839-7846[ISI][Medline]
Langen R., I. J. Chang, J. P. Germanas, J. H. Richards, J. R. Winkler, H. B. Gray, 1995 Electron tunneling in proteinscoupling through a beta strand Science 268:1733-1735[ISI][Medline]
Lazcano A., S. L. Miller, 1996 The origin and early evolution of lifeprebiotic chemistry, the pre-RNA world, and time Cell 85:793-798[ISI][Medline]
Madden T. L., R. L. Tatusov, J. Zhang, 1996 Applications of network BLAST server Meth. Enzymol 266:131-141[ISI][Medline]
Marcus R. A., N. Sutin, 1985 Electron transfers in chemistry and biology Biochim. Biophys. Acta 811:265-322[ISI]
Medvedev D. M., I. Daizadeh, A. A. Stuchebrukhov, 2000 Electron transfer tunneling pathways in bovine heart cytochrome c oxidase J. Am. Chem. Soc 122:6571-6582[ISI]
MSI. 1995 Insight II User Guide, October 1995 Biosym/MSI, San Diego
Mouesca J.-M., J. L. Chen, L. Noodleman, D. Bashford, D. A. Case, 1994 Density functional/Poisson-Boltzmann calculations of redox potentials for iron-sulfur clusters J. Am. Chem. Soc 116:11898-11914[ISI]
Nealson K. H., P. G. Conrad, 1999 Life: past, present and future Philos. Trans. R. Soc. Lond. B 354:1923-1939[ISI][Medline]
Newton M. D., 1991 Quantum chemical probes of electron-transfer kinetics: the nature of donor-acceptor interactions Chem. Rev 91:767-792[ISI]
Page C. C., C. C. Moser, X. Chen, P. L. Dutton, 1999 Natural engineering principles of electron tunneling in biological oxidation-reduction Nature 402:47-52[ISI][Medline]
Quinkal I., P. Kyritsis, I. Kohzuma, S.-C. Im, A. G. Sykes, J.-M. Moulis, 1994 On the role of conserved proline residues in the structure and function of Clostridium pasteurianum 2[4FE-4S] ferredoxin Prot. Eng 7:681-687[Abstract]
Skourtis S. S., D. N. Beratan, 1999 Theories of structure-function relationships for bridge-mediated electron transfer reactions. Electron transfer-from isolated molecules to biomolecules, pt 1 Adv. Chem. Phys 106:377-452[ISI]
Stuchebrukhov A. A., 1996 Tunneling currents in electron transfer reactions in proteins J. Chem. Phys 104:8424-8432[ISI]
Stuchebrukhov A. A., 1996a. Tunneling currents in electron transfer reaction in proteins 2. Calculation of electronic superexchange matrix element and tunneling currents using nonorthogonal basis sets J. Chem. Phys 105:10819-10829.[ISI]
Stuchebrukhov A. A., 1997b. Tunneling currents in proteins: nonorthogonal atomic basis sets and Mulliken population analysis J. Chem. Phys 107:6495-6498[ISI]
Walker J. E., 1992 The NADH-Ubiquinone oxidoreductase (complex-I) of respiratory chains Q. Rev. Biophys 25:253-324[ISI][Medline]
Whangbo M. H., M. Evain, T. Hughbanks, M. Kertesz, S. Wijcycsekera, C. Wilker, C. Zheng, R. Hoffmann, 1987 Extended Hückel molecular and crystal calculations EHMACC was obtained from the Quantum Chemistry Program Exchange (QCPE) number 571
Woese C. R., 1982 Archaebacteria and cellular origins: an overview Zentralbl Bakteriolog. Mikrobiol. Hyg. Ser. C 3:1-17
Yates K., 1978 Hückel molecular orbital theory Academic Press, New York