Sub-grid and small scale processes occur in various ecosystems and landscapes (e.g., periglacial ecosystems, peatlands and vegetation patterns). These local heterogeneities are often important or even fundamental to better understand general and large scale properties of the system, but they are either ignored or poorly parameterized in regional and global models. Because of their small scale, the underlying generating processes can be well explained and resolved only by local mechanistic models, which, on the other hand, fail to consider the regional or global influences of those features. A challenging problem is then how to deal with these interactions across different spatial scales, and how to improve our understanding of the role played by local soil heterogeneities in the climate system. This is of particular interest in the northern peatlands, because of the huge amount of carbon stored in these regions. Land-atmosphere greenhouse gas fluxes vary dramatically within these environments. Therefore, to correctly estimate the fluxes a description of the small scale soil variability is needed. Applications of statistical physics methods could provide useful tools to upscale local features of the landscape, relating them to large-scale properties. To test this approach we considered a case study: the polygonal tundra. Cryogenic polygons, consisting mainly of elevated dry rims and wet low centers, pattern the terrain of many subartic regions and are generated by complex crack-and-growth processes. Methane, carbon dioxide and water vapor fluxes vary largely within the environment, as an effect of the small scale processes that characterize the landscape. It is then essential to consider the local heterogeneous behavior of the system components, such as the water table level inside the polygon wet centers, or the depth at which frozen soil thaws. We developed a stochastic model for this environment using Poisson-Voronoi diagrams, which are able to upscale statistical large scale properties of the system taking into account the main processes within the single polygons. We then compare the results with available recent field studies and demonstrate that the model captures the main statistical characteristics of the landscape and describes its dynamical behavior under climatic forcings (e.g., precipitation and evapotranspiration). In particular, we model and analyze water table dynamics, which directly influences greenhouse gas emissions and changes in the system. Hydraulic interconnectivities and large-scale drainage may also be investigated through percolation properties and thresholds in the Voronoi-Deleaunay graph.