  We present a factorization theorem valid near the kinematic threshold
$z=Q^2/\hat{s}\to 1$ of the partonic Drell-Yan process $q\bar q\to\gamma^*+X$
for general subleading powers in the $(1-z)$ expansion. We then consider the
specific case of next-to-leading power. We discuss the emergence of collinear
functions, which are a key ingredient to factorization starting at
next-to-leading power. We calculate the relevant collinear functions at
$\mathcal{O}(\alpha_s)$ by employing an operator matching equation and we
compare our results to the expansion-by-regions computation up to the
next-to-next-to-leading order, finding agreement. Factorization holds only
before the dimensional regulator is removed, due to a divergent convolution
when the collinear and soft functions are first expanded around $d=4$ before
the convolution is performed. This demonstrates an issue for threshold
resummation beyond the leading-logarithmic accuracy at next-to-leading power.
