Internal gravity waves induce a horizontal mean flow that interacts with the waves if they are of moderately large amplitude. Previous studies of finite-amplitude Boussinesq wavepackets have shown that this interaction results in modulational instability for waves with frequency omega > 2^{-1/2} N. Otherwise, the waves are modulationally stable, and their amplitude decreases faster than predicted by the theory of linear dispersion. In this work, a theoretical form for the wave-induced mean flow of non-Boussinesq internal gravity waves and the corresponding weakly nonlinear Schrodinger equation will be presented. The solution of the Schrodinger equation is numerically integrated for horizontally periodic, vertically localized wavepackets. We consider Gaussian wavepackets with small initial amplitude, propagating upwards through a background in which density decreases exponentially with height. As they propagate, the wave amplitude increases and weakly nonlinear effects develop. The Schrodinger equation results, when compared with fully nonlinear numerical simulations, are found to capture the dominant characteristics of the wave evolution, in the linear and weakly nonlinear regimes. In particular, hydrostatic waves propagate well above the level at which linear theory predicts they should overturn. Implications for gravity wave drag parameterization in the atmosphere are discussed.