We derive stellar masses from SED fitting to rest-frame optical and UV fluxes for 401 star-forming galaxies at z \sim 4 ,~ { } 5 , and 6 from Hubble-WFC3/IR observations of the ERS combined with the deep GOODS-S Spitzer/IRAC data ( and include a previously-published z \sim 7 sample ) . A mass-luminosity relation with strongly luminosity-dependent \mathcal { M } /L _ { UV } ratios is found for the largest sample ( 299 galaxies ) at z \sim 4 . The relation \mathcal { M } \proptoL _ { UV, 1500 } ^ { 1.7 ( \pm 0.2 ) } has a well-determined intrinsic sample variance of 0.5 dex . This relation is also consistent with the more limited samples at z \sim 5 - 7 . This z \sim 4 mass-luminosity relation , and the well-established faint UV luminosity functions at z \sim 4 - 7 , are used to derive galaxy mass functions ( MF ) to masses \mathcal { M } \sim 10 ^ { 8 } at z \sim 4 - 7 . A bootstap approach is used to derive the MFs to account for the large scatter in the \mathcal { M } –L _ { UV } relation and the luminosity function uncertainties , along with an analytical crosscheck . The MFs are also corrected for the effects of incompleteness . The incompleteness-corrected MFs are steeper than previously found , with slopes \alpha _ { M } \sim - 1.4 to -1.6 at low masses . These slopes are , however , still substantially flatter than the MFs obtained from recent hydrodynamical simulations . We use these MFs to estimate the stellar mass density ( SMD ) of the universe to a fixed M _ { UV,AB } < -18 as a function of redshift and find a SMD growth \propto ( 1 + z ) ^ { -3.4 \pm 0.8 } from z \sim 7 to z \sim 4 . We also derive the SMD from the completeness-corrected MFs to a mass limit \mathcal { M } \sim 10 ^ { 8 } M _ { \odot } . Such completeness-corrected MFs and the derived SMDs will be particularly important for model comparisons as future MFs reach to lower masses .