We present a formulation of Raman spectroscopy in molecular junctions based on a many-body state representation of the molecule. The approach goes beyond the previous effective single orbital formalism and provides a convenient way to incorporate computational methods and tools proven for equilibrium molecular spectroscopy into the realm of current carrying junctions. The presented framework is illustrated by first principle simulations of Raman response in a three-ring oligophenylene vinylene terminating in amine functional groups (OPV3) junction. The calculated shift in Stokes lines and estimate of vibrational heating by electric current agree with available experimental data. In particular, our results suggest that participation of the OPV3 cation in Raman scattering under bias may be responsible for the observed shift, and that the direction of the shift depends on renormalization of normal modes. This work is a step toward atomistic quantum ab initio modeling of the optical response of nonequilibrium electronic dynamics in molecular junctions.