  We describe the perturbative calculation of the transverse parton
distribution functions in all partonic channels up to next-to-next-to-leading
order based on a gauge invariant operator definition. We demonstrate the
cancellation of light-cone divergences and show that universal
process-independent transverse parton distribution functions can be obtained
through a refactorization. Our results serve as the first explicit higher-order
calculation of these functions starting from first principles, and can be used
to perform next-to-next-to-next-to-leading logarithmic $q_T$ resummation for a
large class of processes at hadron colliders.
