We present numerical scaling results for the energy level statistics in orthogonal and symplectic tight-binding Hamiltonian random matrix ensembles defined on disordered two and three-dimensional electronic systems with and without spinorbit coupling (SOC), respectively. In the metallic phase for weak disorder the nearest level spacing distribution functionP(S), the number variance , and the two-point correlation functionK2(e), are shown to be described by the Gaussian random matrix theories. In the insulating phase, for strong disorder, the correlations vanish for large scales and the ordinary Poisson statistics is asymptotically recovered, which is consistent with localization of the corrosponding eigenstates. At the Anderson metal-insulator transition we obtain new universal scale-invariant distribution functions describing the critical spectral density fluctuations.