Numerical integration of the double- and triple-box integrals using threshold subtraction

We numerically integrate finite two- and three-loop scalar integrals using the threshold subtraction method. This represents a first step towards extending our calculation of the $N_f$-part to the full NNLO virtual corrections for the production of massive vector bosons.