  HI intensity mapping is a new observational technique to map fluctuations in
the large-scale structure of matter using the 21 cm emission line of atomic
hydrogen (HI). Sensitive radio surveys have the potential to detect Baryon
Acoustic Oscillations (BAO) at low redshifts (z < 1) in order to constrain the
properties of dark energy. Observations of the HI signal will be contaminated
by instrumental noise and, more significantly, by astrophysical foregrounds,
such as Galactic synchrotron emission, which is at least four orders of
magnitude brighter than the HI signal. Foreground cleaning is recognised as one
of the key challenges for future radio astronomy surveys. We study the ability
of the Generalized Needlet Internal Linear Combination (GNILC) method to
subtract radio foregrounds and to recover the cosmological HI signal for a
general HI intensity mapping experiment. The GNILC method is a new technique
that uses both frequency and spatial information to separate the components of
the observed data. Our results show that the method is robust to the complexity
of the foregrounds. For simulated radio observations including HI emission,
Galactic synchrotron, Galactic free-free, radio sources and 0.05 mK thermal
noise, we find that we can reconstruct the HI power spectrum for multipoles 30
< l < 150 with 6% accuracy on 50% of the sky for a redshift z ~ 0.25.
