Land surface interacts strongly with the atmosphere at all scales. This has a considerable impact on the hydrologic cycle and the climate. Therefore, in order to produce realistic simulations with climate models, their land-surface processes must be parameterized accurately. Because continental surfaces are usually extremely heterogeneous over the resolvable scales considered in these models, surface parameterizations based on the 'big leaf-big stoma' approach (that assume grid-scale homogeneity) fail to represent the land-atmosphere interactions that occur at much smaller scales. A parameterization based on a statistical-dynamical approach is suggested here. With this approach, each surface grid element of the numerical model is divided into homogeneous land patches (i.e., patches with similar internal heterogeneity). Assuming that horizontal fluxes between the different patches within a grid element are small as compared to the vertical fluxes, patches of the same type located at different places in the grid can be regrouped into one subgrid surface class. Then, for each one of the subgrid surface classes, probability density functions (pdf) are used to characterize the variability of the different parameters of the soil-plant-atmosphere system. These pdf are combined with the equations of the model that describe the dynamic and the energy and mass conservations in the atmosphere. The potential application of this statistical-dynamical parameterization is illustrated by simulating (i) the development of an agricultural area in an arid region and (ii) the process of deforestation in a tropical region. Both cases emphasize the importance of land-atmosphere interactions on regional hydrologic processes and climate.