Numerical Renormalization Group for Bosonic Systems and Application to the Subohmic Spin-Boson Model