  We present an \textit{ab initio} theory for superconductors, based on a
unique mapping between the statistical density operator at equilibrium, on the
one hand, and the corresponding one-body reduced density matrix $\gamma$ and
the anomalous density $\chi$, on the other. This new formalism for
superconductivity yields the existence of a universal functional
$\mathfrak{F}_\beta[\gamma,\chi]$ for the superconductor ground state, whose
unique properties we derive. We then prove the existence of a Kohn-Sham system
at finite temperature and derive the corresponding Bogoliubov-de Gennes-like
single particle equations. By adapting the decoupling approximation from
density functional theory for superconductors we bring these equations into a
computationally feasible form. Finally, we use the existence of the Kohn-Sham
system to extend the Sham-Schl\"uter connection and derive a first
exchange-correlation functional for our theory. This reduced density matrix
functional theory for superconductors has the potential of overcoming some of
the shortcomings and fundamental limitations of density functional theory of
superconductivity.
