An inverse dynamic procedure for spatial multibody systems containing flexible bodies is developed in the relative joint coordinate space. Constraint acceleration equations are derived in terms of relative coordinates using the velocity transformation technique. An inverse velocity transformation operator, which transforms the Cartesian velocities to the relative velocities, is derived systematically corresponding to the types of kinematic joints connecting the bodies and the system reference matrix. Using the resulting matrix, the joint reaction forces and moments are analyzed in the Cartesian coordinate space. The formulation is illustrated by means of two numerical examples.