  Given a pair of self-adjoint operators $H$ and $V$ such that $V$ is bounded
and $(H+V-i)^{-1}-(H-i)^{-1}$ belongs to the Schatten-von Neumann ideal
$\mathcal{S}^n$, $n\ge 2$, of operators on a separable Hilbert space, we
establish higher order trace formulas for a broad set of functions $f$
containing several major classes of test functions and also establish existence
of the respective locally integrable real-valued spectral shift functions
determined uniquely up to a low degree polynomial summand. Our result
generalizes the result of \cite{PSS13} for Schatten-von Neumman perturbations
$V$ and settles earlier attempts to encompass general perturbations with
Schatten-von Neumman difference of resolvents, which led to more complicated
trace formulas for more restrictive sets of functions $f$ and to analogs of
spectral shift functions lacking real-valuedness and/or expected degree of
uniqueness. Our proof builds on a general change of variables method derived in
this paper and significantly refining those appearing in \cite{vNS21,PSS15,S17}
with respect to several parameters at once.
