  We derive the statistics of the time-delay matrix (energy derivative of the
scattering matrix) in an ensemble of superconducting quantum dots with chaotic
scattering (Andreev billiards), coupled ballistically to $M$ conducting modes
(electron-hole modes in a normal metal or Majorana edge modes in a
superconductor). As a first application we calculate the density of states
$\rho_0$ at the Fermi level. The ensemble average
$\langle\rho_0\rangle=\delta_0^{-1}M[\max(0,M+2\alpha/\beta)]^{-1}$ deviates
from the bulk value $1/\delta_0$ by an amount depending on the
Altland-Zirnbauer symmetry indices $\alpha,\beta$. The divergent average for
$M=1,2$ in symmetry class D ($\alpha=-1$, $\beta=1$) originates from the
mid-gap spectral peak of a closed quantum dot, but now no longer depends on the
presence or absence of a Majorana zero-mode. As a second application we
calculate the probability distribution of the thermopower, contrasting the
difference for paired and unpaired Majorana edge modes.
