An integrated computational methodology for interpreting NMR spin relaxation in proteins has been developed. It combines a two-body coupled-rotator stochastic model with a hydrodynamics-based approach for protein diffusion, together with molecular dynamics based calculations for the evaluation of the coupling potential of mean force. The method is applied to ¹⁵N relaxation of N-H bonds in the Rho GTPase binding (RBD) domain of plexin-B1, which exhibits intricate internal mobility. Bond vector dynamics are characterized by a rhombic local ordering tensor, S, with principal values S₀² and S₂², and an axial local diffusion tensor, D₂, with principal values D(2,||) and D(2,⊥). For α-helices and β-sheets we find that S₀² ~ -0.5 (strong local ordering), -1.2 < S₂² < -0.8 (large S tensor anisotropy), D(2,⊥) ~ D₁ = 1.93 × 10⁷ s⁻¹ (D₁ is the global diffusion rate), and log(D(2,||)/D₁) ~ 4. For α-helices the z-axis of the local ordering frame is parallel to the C(α)-C(α) axis. For β-sheets the z-axes of the S and D₂ tensors are parallel to the N-H bond. For loops and terminal chain segments the local ordering is generally weaker and more isotropic. On average, D(2,⊥) ~ D₁ also, but log(D(2,||)/D₁) is on the order of 1-2. The tensor orientations are diversified. This study sets forth an integrated computational approach for treating NMR relaxation in proteins by combining stochastic modeling and molecular dynamics. The approach developed provides new insights by its application to a protein that experiences complex dynamics.