In this thesis, we observationally parametrize the average properties of the stellar haloes of galaxies in order to constrain the physics of galaxy formation. By stacking aligned mosaics of a large number of face-on central galaxies from the Sloan Digital Sky Survey (SDSS), we study the properties of the stellar haloes of galaxies as a function of stellar mass (from MPA-JHU catalogue), halo mass (from Yang et al. 2007 group catalogue) and morphology. Using stacks of mock images of galaxies from the Illustris simulations, we show that the outer light fraction derived from fitting double S\'{e}rsic models to the 2-D surface brightness distribution of galaxy stacks provides an upper limit (within 0.1 dex) of the mean accreted stellar mass fraction. For the SDSS stacks, we find that the outer light fraction is a function of stellar mass and galaxy type, increasing from 30% to 70% and from 2% to 25% for early and late type galaxies respectively over the mass range between $10^{10.0} \msun$ to $10^{11.4} \msun$. Above the characteristic mass ($\log M_{halo} \sim 12.5$), we find that the outer light fraction is a stronger function of halo mass than stellar mass. Below the characteristic mass, the outer light fraction is a strong function of a galaxy concentration ($R_{90}/R_{50}$). We further parametrize the surface brightness distribution of the stellar halo of the galaxy stack by estimating its outer slope and ellipticity. We find that the outer slope is an increasing function of halo mass and concentration, while the outer ellipticity is a increasing function of stellar mass and concentration. The $g$-$r$ colour of the stellar population in the stellar halo is bluer than in the main galaxy, and the colour of the stellar halo is redder for higher mass galaxies. We find that our observational constraints agree well with the Illustris simulations above the characteristic mass. However, the simulations fail to reproduce the data below the characteristic mass where the accreted mass fraction is an increasing function of stellar mass. Using our SDSS stacks, we derive average flux corrections to the SDSS \texttt{Model} magnitudes, finding corrections ranging from 0.05 to 0.32 mag for the highest stellar mass galaxies. We apply these corrections to the MPA-JHU stellar masses for a complete sample of half a million galaxies from the SDSS survey to derive a corrected galaxy stellar mass function at $z=0.1$ in the stellar mass range $9.5<\log(M_\ast/M_\odot)<12.0$. We find that the flux corrections and the use of the MPA-JHU stellar masses have a significant impact on the massive end of the stellar mass function, making the slope significantly shallower than that estimated by Li \& White (2009), but steeper than derived by Bernardi et al. (2013). This corresponds to a mean comoving stellar mass density of galaxies with stellar masses $\log(M_\ast/M_\odot) \ge 11.0$ that is a factor of 3.36 larger than the estimate by Li \& White (2009), but is 43\% smaller than reported by Bernardi et al. (2013).