In the present work, density functional theory (DFT) has been used to investigate CO binding to the hexacoordinated heme in neuroglobin (Ngb) protein. Structural relaxation of the selected model system in the protein environment has been fully included by the alternative quantum and molecular mechanical optimizations. The polarized continuum model (PCM) was used to simulate interaction between the model system and the protein environment. The CO binding could take place in a concerted way and a barrier of 17.9 kcal mol(-1) was predicted on the concerted singlet pathway, which is not favorable in energy. The adiabatically sequential pathway requires an energy of 14.5 kcal mol(-1) for formation of the singlet intermediate. There exist two nonadiabatic sequential pathways for the CO binding, which involves the triplet and quintet states of intermediate. Both the singlet/triplet and singlet/quintet intersections play an important role in nonadiabatic sequential processes, which enhance the probability that the processes occur. The nonadiabatic processes that involve the triplet and quintet states of intermediate are the most probable pathways for the CO binding to the hexacoordinated heme in Ngb to form the product complex.