This work proposes a methodology for in situ parameter identification using system-level measurements of (flexible) multibody systems, opposed to dedicated component-level identification. The sensitivity information employed for the optimization is obtained using the adjoint variable method (AVM). This method has the advantage of obtaining sensitivity information at a computational cost independent of the amount of model parameters. The underlying flexible multibody formulation employed is a novel approach called the flexible natural coordinates formulation (FNCF). This formulation combines the advantageous properties of the floating frame of reference formulation (FFRF) and the generalized component mode synthesis (GCMS) methods and results in a constant mass and stiffness matrix with quadratic constraint equations. This work shows how the specific structure of equations obtained through FNCF drastically reduces the complexity of the AVM as the simulation derivatives can be readily obtained and are of limited order. The proposed approach has been implemented in an in-house object-oriented matlab multibody code. The methodology is illustrated by identifying 13 model parameters of a MacPherson suspension model, in situ and using system-level measurements.