The reasons that the temperature dependence of the NMR isotropic shifts of model ferrihemes and ferriheme proteins deviate from Curie behavior have been analyzed by considering the energies of the valence orbitals of the metal and the porphyrinate. For low-spin Fe(III), overlap of the e-symmetry pi orbitals of a symmetrical porphyrin ring and the d(pi) orbitals of the metal produces two low-energy molecular orbitals that are mainly porphyrin in character and are filled and two high-energy (valence) molecular orbitals that are mainly metal in character and contain three electrons. The odd electron in the valence set thus gives rise to the spin delocalization that results in the observed contact shift pattern of these systems. Unsymmetrical substitution and/or presence of a planar axial ligand that is prevented from rotation removes the degeneracy of these e(pi) orbitals, producing a system in which the energy separation between the two formerly degenerate pi orbitals, Delta E(pi), is of the order of only tens of cm(-1) for the former or quite large (several times k(B)T) for the latter. In either case, both orbitals are utilized for spin delocalization to a significant extent as the temperature is varied, according to their varying Boltzmann populations. Such a two-level system obeys a modified Curie law that takes into account the thermal population of the two levels as a function of temperature. In fact, the temperature dependence of some of the contact shifts of model hemes or heme proteins may show anti-Curie behavior if Delta E(pi) is large compared to k(B)T at ambient temperatures. Such anti-Curie behavior has been observed for two of the heme methyl resonances of several cytochromes c and b(5) and cyanometmyoglobins or -hemoglobins, where the axial methionine pi-symmetry lone pair or histidine imidazole plane orientation, respectively, is believed to be the important factor in determining Delta E(pi). Assuming reasonable energy separations of the two valence e(pi) orbitals, from very small to quite large (similar to 1000 cm(-1)), the expected temperature dependence of the contact shifts has been calculated for an assumed set of valence MO coefficients. These results have then been compared to the experimental isotropic shifts of several model heme systems having unsymmetrical substitution patterns and/or one fixed axial ligand and to several heme proteins. Using a computer program developed to fit the observed isotropic shifts to the two-level equation, Delta E(pi) was estimated from the temperature dependence of the isotropic shifts of the protons of the beta-pyrrole substituents of the above-mentioned systems. In the case of the proteins investigated, Aplysia cyanometmyoglobin and cytochrome b(5), the values of Delta E(pi), obtained from analysis of proton isotropic shifts are similar to those calculated from EPR g values measured at low temperatures, while for model hemins, the values of Delta E(pi) obtained are smaller than those predicted and vary in accord with the expectations as to the rigidity, or lack thereof, of the orientation of at least one planar axial ligand, indicating that thermal averaging of the two levels due to rapid rotation (or libration) of the axial ligand is fast on the NMR time scale. This same two-level approach could be applied to any system in which there is a thermal equilibrium between two states separated by an energy within several factors of k(B)T.