The free energies of two low pressure ice forms are calculated over a wide range of temperatures in order to explain relative stability and their negative thermal expansivities in the low temperature regime. One-hundred proton-disordered configurations for hexagonal and cubic ices are generated. The Helmholtz free energy is approximated to a sum of the minimum potential energy, the harmonic free energy and the configurational entropy arising from the disordered nature of protons. The Gibbs free energy at a given temperature is minimized with respect to the volume of the system. This enables us to evaluate the thermal expansivity at a fixed temperature and pressure from only intermolecular interaction potentials. The negative thermal expansivity of ices in low temperature is successfully reproduced. This arises mainly from the bending motion of hydrogen bonded molecules.