Evaluation of Snow Depth and Snow Cover Fraction Simulated by Two Versions of the Flexible Global Ocean-Atmosphere-Land System Model

The land surface area covered by snow varies greatly through the year over the Northern Hemisphere (NH). This seasonal variation in snowpack is treated as a significant feature of the climate because snow has special properties such as high reflectivity and low thermal conductivity (Essery, 1997), which can influence the exchange of energy and water between atmosphere and land and eventually lead to changes in atmosphere circulation and climate. Knowledge of the relationships between snow cover and climate can be dated back to the 19th century (Voeikov, 1889). From then on, numerous statistical and modeling studies have been carried out over regional/continental scales. Although not exactly the same conclusions are always obtained, all research in this area emphasizes the significant impact of snow cover on climate (Cohen and Rind, 1991; Popova, 2007; Peings and Douville, 2010). In view of the relationship between accumulated Eurasian snow during winter and spring and East Asian rainfall the following summer, Eurasian snow is considered as one of the main predictors for precipitation during the rainyseason (Wu and Kirtman, 2007). Therefore, snow studies of this kind are extremely important on a number of levels.

At present, general circulation models (GCMs) are accepted as one of the most effective tools in helping us to understand the past and present climate, as well as predict its future changes (Shackley et al., 1998). Snow has been included in GCMs due to its importance in climate change, and accurate simulation of snow has a significant impact on the prediction of weather and climate. Moreover, snow is a particularly good diagnostic element for model evaluation, since correct representations of snow cover and snow depth rely closely upon accurate simulations of both temperature and precipitation. As most GCMs have used snow schemes in their land surface models (LSMs), some specific snow comparison projects have been carried out:

The most noticeable of these projects was phase 2(d) of the Project for Intercomparison of Land Surface Parameterization Schemes (PILPS), which showed that land surface schemes (LSSs), as a group, can capture the general patterns of accumulation and ablation on an interannual basis, but systematic scatter exists during the early part of the snow season and ablation events (Slater et al., 2001). Other projects include the first and second phases of the Snow Models Intercomparison Project (SnowMIP), in which 26 and 33 land surface models from countries including Russia, France and China were respectively included under the same atmospheric forcing conditions. All snow schemes included in PILPS 2(d) and SnowMIP were validated in offline simulations over several site locations. The objectives of intercomparison projects do not involve choosing one best snow model, but instead aim to compare the divergence of various snow schemes by considering the same atmospheric forcing conditions as the input, and further identify the key processes for each application. It has been found that most models compared in this way perform well in terms of simulating snow accumulation, but vary greatly in their representations of snowmelt. Furthermore, most models show faster accumulation compared to observations (Etchevers et al., 2002, 2004; Rutter et al., 2009). However, in addition to the differences in simulated snow caused by the various snow schemes, interactions between atmosphere and land also play important roles in the simulation of snow by GCMs. This aspect has not been considered by the intercomparison projects and studies carried out to date.

The two most recent versions of the Flexible Global Ocean-Atmosphere-Land System (FGOALS) models are the Grid-point Version 2 (g2) (Li et al., 2013b) and the Spectral Version 2 (s2) (Bao et al., 2013). Both participated in the fifth Coupled Model Intercomparison Project (CMIP5), and both include exactly the same model components for the ocean, sea ice and land, but differ in their atmospheric components. This provides a very good opportunity to study the influence of the atmospheric component on the simulation of snow. Therefore, the main purpose of the present reported study was to evaluate the performance of the two models in simulating snow cover, and try to analyze the reasons leading to any differences found.

FGOALS-g2 and FGOALS-s2 used in this study are fully coupled climate models. Both incorporate four components, including the same ocean circulation model (LICOM2) (Liu et al., 2012; Lin et al., 2013), the same land surface process model (CLM3) (Oleson et al., 2004), and the same sea ice model (CICE4-LASG). However, the atmospheric circulation models are different: for FGOALS-g2 it is Grid Atmospheric Model of IAP/LASG, Version 2 (GAMIL2), but for FGOALS-s2 it is Spectral Atmospheric Circulation Model of IAP/LASG, Version 2 (SAMIL2). All components are coupled together by the CPL6 coupler from the National Center for Atmospheric Research (NCAR) (Craig et al., 2005). For convenience, FGOALS-g2 and FGOALS-s2 are hereafter abbreviated as g2 and s2, respectively. Both models have been used to perform CMIP5 experiments, and the results are accessible to researchers throughout the world via the State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics (LASG) node of the Earth System Grid (ESG).

Both GAMIL2 and SAMIL2 have the same dynamic core (Wang et al., 2004). GAMIL2 employs a hybrid horizontal grid with a Gaussian grid of 2.8^° between 65.58^°S and 65.58^°N and a weighted equal-area grid poleward of 65.58^°. SAMIL2 is a spectral transform model with a horizontal resolution of R42, which is a resolution of approximately 1.66^° (lat) ×2.81^° (lon). The main physical difference between the two atmospheric models is their cloud-related processes. Compared to s2, g2 has a two-moment cloud microphysics scheme adopting two different physical-based aerosol activation parameterizations is employed to estimate the aerosol indirect effects (Rasch and Kristjansson, 1998; Morrison and Gettelman, 2008; Shi et al., 2010). A diagnostic Slingo-type scheme (Slingo, 1987) is employed to calculate the cloud fraction that depends on relative humidity, water vapor, atmospheric stability and convective mass fluxes. In s2, the radiation is calculated by adopting the modified Edwards-Slingo scheme (Edwards and Slingo, 1996; Sun and Rikus, 1999a, 1999b), which is able to accommodate the direct effects of aerosols, but cannot represent indirect aerosol effects. The cloud is diagnosed by a scheme based on vertical motion and relative humidity (Liu and Wu, 1997). Detailed information about GAMIL2 and SAMIL2 can be found in (Li et al., 2013a) and (Bao et al., 2013).

The land surface model adopted by g2 and s2 is CLM3.0, which shares the same grid as the atmospheric model. CLM3.0 is constructed with 10 uniformly-distributed vertical soil layers, with (at most) five snow layers, and one vegetation layer. The flux is calculated separately for different underlying surfaces by considering the inhomogeneity on the sub-grid scale of land-surface characteristics. The snow model in CLM3 is a multilayer physical-based model, primarily simplified from (Anderson, 1976) and (Jordan, 1991) by parameterizing gravitational liquid water flow and grain size growth while neglecting the water vapor phase. The state variables for snow are partial volume of liquid water and ice, snow density, and temperature. Total snow mass is divided or combined at every time step in response to changes in layer depth because of snowfall, sublimation and ablation while conserving energy and mass. Three mechanisms for changing snow characteristics are included: destruction, overburdening, and melting. For detailed information, the reader is referred to (Oleson et al., 2004).

LICOM is the oceanic component of g2 and s2, which is the fourth generation of the LASG's Ocean General Circulation Model (OGCM) (Liu et al., 2004a, 2004b), the first generation having been developed at (LASG/ Institute of Atmospheric Physics (IAP) in 1989 (Zhang and Liang, 1989). The model's horizontal resolution is 1^°× 1^° for most areas, but within the 2^° around the Equator in the meridional direction it is 0.5^°× 0.5^°; its vertical resolution has been increased to 10 m for the uppermost 150 m compared with the initial version. Details concerning the improvements of the physical processes of LICOM are covered in (Lin et al., 2007) and (Liu et al., 2012). The Los Alamos sea ice model CICE4-LASG, which shares the same grid system as that of the ocean component (LICOM2), is an updated version of Community Sea Ice Model (CSIM) version 4 (Collins et al., 2006). Detailed information about the design of CICE4-LASG and its performance in CMIP5 are covered in detail in (Xu et al., 2012).

The monthly mean outputs from the historical simulations of the CMIP5 experiment were the focus of this study. Four (g2) and three (s2) groups of historical simulations were performed as ensemble members to reduce the uncertainties arising from the different initial conditions. The ensemble-based historical experiments were carried out based on different beginning years in the pre-industrial control experiments. In our discussion of the results, below, we focus mainly on comparing the two models based on their simulations over the NH for the period 1956-2005, during which more observational data are available. In addition, 50 years is long enough to determine the year-to-year climate variability of land variables according to (Manabe and Stouffer, 1996).

Observational datasets of snow depth (SD) and snow cover fraction (SCF) were used to validate the models. Observational elements of precipitation and surface temperature were also used.

3.1.1. Snow depth

The global snow depth climatology of the U.S. Air Force/Environmental technical Application Center (USAF/ ETAC) (Foster and Davy, 1988) is considered to be one of the most credible sources of SD data, and was therefore employed in the present study for validation against modeled SDs. The monthly USAF/ETAC global snow depth climatology is obtained based on a set of station-referenced measurements of snow depth and has been used in several studies to validate snow models (Douville et al., 1995a, 1995b; Niu and Yang, 2006).

3.1.2. Snow cover fraction

Both the National Oceanic and Atmospheric Administration (NOAA) visible data and the Moderate-resolution Imaging Spectroradiometer (MODIS) monthly SCF product were used in the present study. Monthly SCF of the NH was derived from weekly values of the "weekly digital NH snow and ice product", which has been produced by the NOAA and the National Environmental Satellite Data and Information Service (NESDIS) since 1966. Generally, the NOAA charts are recognized as the most suitable method to obtain snow cover information on the global scale. In addition, an initial version of the gridded global monthly snow cover product from 1971 to 1994 was also used in this study.

3.1.3. Temperature

The Willmott-Matsuura climatology of monthly terrestrial air temperature (version 3.01) (Willmott and Matsuura, 2000) was used, the data for which were obtained from the Legates and Willmott (1990a, 1990b) station archives of monthly air temperature and the Global Historical Climatology Network (GHCN). Station data were interpolated to a grid field with a spatial resolution of 0.5^°× 0.5^° by using a distance weighted method and temperature was adjusted by considering the effect of elevation (Willmott and Robeson, 1995). Furthermore, CRUTEMP3 was applied to validate the long-term temperature variations from 1850 to 2005. The CRUTEMP3 dataset is based on records of monthly average temperature obtained from more than 4000 weather stations; temperature anomalies were calculated for every station from the 1961-90 average.

3.1.4. Precipitation

The CPC Merged Analysis of Precipitation (CMAP) dataset, which is compiled from a variety of satellite and rain gauge data, was employed. A detailed description of this dataset can be found in Xie and Arkin (1996, 1997).

All the observational data and simulated results used in this work were re-gridded to 128 (lon) ×60 (lat) using the bilinear interpolation method. Basic statistical methods and error analysis were also employed to evaluate the SD and SCF simulated by the two models.

3.2.1. Spatial correlation

The correlation coefficient is a statistical index able to reflect the closeness of the relationship between different variables. The spatial correlation coefficient is defined as the covariance divided by the standard deviations of the two variables and is usually used to quantify the degree of consistency between observed and simulated spatial patterns.

3.2.2. Relative standard deviation

Relative standard deviation (RSD) is obtained by dividing the standard deviation by its average, which is also considered as the absolute value of the coefficient of variation (CV) and calculated as a percentage. It can be used in studies like the present one to quantify the interannual variability of every month. However, it is important to be aware that, when the mean value is close to zero, the CV is sensitive to small changes in the mean, then limiting its usefulness.

Figure 1.  Spatial distribution of (a, c, e) mean SD (units: m) and (b, d, f) SCF in winter during 1956-2005: (a, b) g2; (c, d) s2; (e, f) observations.

3.2.3. Root-mean-square error and linear regression analysis

Root-mean-square error (RMSE) is frequently used as a measurement of the differences between values predicted by a model and the observation, and is recognized as a good measurement of the accuracy of simulated results. In the present study, the linear regression coefficient is used to reveal the temporal trends of variables at the decadal scale.

4.1.1. Spatial correlation analysis

Figure 1(a, c, e) shows the spatial distribution of SD from USAF/ETAC and the two models in winter during 1956-2005. Generally speaking, both models capture the SD pattern and the maximum zone in central Siberia is depicted well. However, a big difference exists over the Tibetan Plateau. The simulated SD by g2 over the Tibetan Plateau is close to 1.0 m, while s2 and the observed SD are about 0.2 m. Due to the lack of observations at high altitudes and the inaccurate description of the cryosphere in the model, both the observed and simulated SDs over the Tibetan Plateau show great discrepancies. The monthly variations of spatial correlation coefficients are shown in Fig. 2a, and it is clear that both g2 and s2 have high positive correlation. The average coefficient is about 0.7, and the difference in correlation coefficients between the two models is negligible. Therefore, the spatial pattern of SD can be successfully captured by both models.

4.1.2. Seasonal cycle

Figure 3 shows the ensemble mean climatological annual cycle of SD for eight representative mid-latitude and high-latitude areas for both models and the USAF/ETAC observational data. As can be seen, g2 substantially overestimates SD during the non-summer seasons in all regions apart from eastern Canada, as compared to the observation and s2. These positive deviations are remarkable in late winter and spring, and more specifically, the largest positive anomalies (of up to 0.3 m) are found in May in Northern Europe. Excessive snow accumulation during these periods results in a delay of the simulated dates by which the snow has completely melted. For example, in Northern Europe, the SD given by g2 reaches its minimum in July, while the observational SD falls to its lowest point in June (Fig. 3e), almost one month earlier. In addition, there is a large discrepancy between the results obtained by the two models and observations over the Tibetan Plateau (Fig.3d). However, this is caused by the lack of snow measurements available for high altitude areas such as the plateau region.

Figure 2.  Monthly variations of spatial correlation coefficients between (a) SD and (b) SCF from g2 and s2 and observations.

Figure 3.  Simulated and observed climatological annual cycle of SD for eight selected regions: (a) China (26^°-55^°N, 100^°-125^°E); (b) Europe (37^°-57^°N, 0^°-32^°E); (c) Central Canada (46^°-61^°N, 123^°-97^°W); (d) Tibetan Plateau (30^°-40^°N, 80^°-100^°E); (e) Northern Europe (60^°-70^°N, 5^°-45^°E); (f) eastern Siberia (50^°-66.5^°N, 90^°-140^°E); (g) Alaska (56^°-75^°N, 167^°-141^°W); and (h) eastern Canada (50^°-60^°N, 80^°-55^°W). The solid black line is the observation, the dashed red line is the simulation by g2, and the dotted blue line is the simulation by s2.

Figure 4.  Frequency distribution of SD for g2, s2 and the USAF/ETAC climatology in China (a, b) and eastern Siberia (c, d) during DJF (a, c) and MAM (b, d).

Figure 5.  Seasonal cycles of SCF for (a) China (excluding the Tibetan Plateau), (b) Europe, (c) Central Canada, (d) Tibetan Plateau, (e) Northern Europe, and (f) eastern Siberia. The dashed line represents MODIS data, and the solid line NOAA data.

Figure 6.  Interannual variation of SCF for (a) China (26^°-55^°N, 100^°-125^°E) (excluding the Tibetan Plateau), (b) Northern Europe (60^°-70^°N, 5^°-45^°E), (c) eastern Siberia (50^°-66.5^°N, 90^°-140^°E), and (d) Central Canada (46^°-61^°N, 123^°-97^°W).

Figure 7.  Decadal trends of SCF [units: % (10 yr)^-1] between 0^° and 90^°N in winter during the period 1971-94: (a) g2; (b) s2; and (c) NOAA.

Figure .  Variation of (a) SD (m), (b) precipitation (mm d^-1), and (c) temperature (^°C) in winter during the period 1956-2005 over the mid-high latitudes of the NH. The dotted line is the trend line. The blue and red lines indicate the simulations of g2 and s2, respectively.

RMSE was used to quantify the errors between simulations and observations. The results show that s2 has a smaller RMSE (0.257) than g2 (0.295) over the NH, indicating that s2 produces less error and a better performance in its simulation of SD over the NH than g2.

4.1.3. Frequency distribution

The frequency distributions of simulated and observed SD for winter (December-January-February; DJF) and spring (March-April-May; MAM) are shown in Fig.4. Here, only China and eastern Siberia are chosen for analysis because the results are much the same over all areas, apart from the Tibetan Plateau. Figure 4b reveals that in spring the frequency distribution of SD in China is captured well by s2. In winter, SD has high frequency for those exceeding 0.1 m and low frequency for those less than 0.1 m. However, SD greater than 0.1 m over China obtained from g2 appears more frequently in both spring and winter. According to USAF/ETAC data, the percentages of SD below 0.1 m in winter and spring for all grid boxes are 82.1% and 93.7%, respectively. In other words, in either season, the observed SD less than 0.1 m in China appears much more frequently than in any other areas. The corresponding values for g2 in winter and spring are 59.1% and 70.0%, respectively, while the percentages for s2 are 79.2% and 95.8%, respectively. Clearly, s2 provides results that are much closer in performance to the observation, but the overall results are in good agreement with the overestimation of mean SD over China during winter and spring by both models.

Over eastern Siberia, the maximum SD frequency lies between 0.2 m and 0.4 m in winter and spring. For the USAF climatology, 72.2% of grid cells show SD as being greater than 0.3 m in winter, while for g2 and s2 the percentages are 74.7% and 69.5%, respectively. During spring, excessive SDs are simulated by the two models, with 73.6% of grid cells showing SDs of greater than 0.3 m in g2, while the corresponding values are 22.2% for USAF and 46.8% for s2. The frequency distributions of SD are captured well by the two models in winter, but are overestimated in spring. However, overall-as with China-the results also agree with the overestimation of SD over eastern Siberia (Fig. 3f).

4.2.1. Spatial correlation analysis

Generally speaking, the spatial distribution of mean SCF, which is larger at higher latitudes than at lower latitudes, and larger on the Tibetan Plateau than other areas located at similar latitudes, is well depicted by the two models, as shown in Fig. 1(b, d, f). However, the simulated SCF of both g2 and s2 is less than observed for almost all mid-high latitude areas of the NH. Although the spatial correlation coefficients between the observation and s2 are slightly less than those with g2, the monthly average of spatial correlation coefficients between the two models and observations is close to 0.8 (Fig 2b), which is larger than that for SD. Therefore, both models are able to simulate the spatial pattern of SCF well.

4.2.2. Seasonal cycle

Although there are uncertainties due to the incomplete knowledge of atmospheric states and the difficulty in distinguishing snow from clouds and so on, agreement between various remotely-sensed data suggests that satellite products comprise one of the most appropriate tools for model validation (Roesch, 2006). It has been found that National Oceanic and Atmospheric Administration (NOAA) and MODIS SCF data are in good agreement with each other over areas of latitude lower than 70^°N, and thus we only evaluated SCF simulations in that region.

Generally speaking, the SCF modeled by g2 is always higher than that given by s2 over the non-summer seasons, especially during late winter and early spring (Fig. 5). This means that the SCF simulated by g2 has more significant seasonal variation compared to s2. Compared with SCF observations from both NOAA and MODIS, s2 underestimates SCF over almost all of the six areas, while g2 produces under- or overestimations of SCF in different regions. For instance, over China (excluding the Tibetan Plateau) (Fig. 5a), g2 overestimates SCF during most of the snow seasons, but over areas of Europe, Central Canada, and Northern Europe (Figs. 5b, c and e), good agreement can be found. For higher altitude areas such as the Tibetan Plateau (Fig. 5d), due to the complex terrain and paucity of observation stations, not much is known yet about the physical processes of these locales. Indeed, all results show great uncertainties for the Tibetan Plateau region. Also, results show that g2 has a smaller RMSE (0.182) compared to s2 (0.211) over the NH, indicating that g2 produces less error and performs better in its simulations of SCF over the NH compared to s2.

4.2.3. Interannual variability

We also examined the models' performance in terms of the interannual variability of SCF, and RSD was used for this by analyzing it for every month. Figure 6 shows the interannual variability of monthly SCF as simulated by both g2 and s2. The regions here are identical to those in Figs. 5a, c, e and f. As all four regions belong to seasonal snow areas, where snow usually melts completely during summer and late spring, and accumulates slowly during early autumn, particular attention is paid to the seasons of winter, late autumn and early spring.

It can clearly be seen that small RSD usually appears in winter, while larger RSD appears in May and October, over all four regions. The inference of this result is that the snowpack is quite thick and stable in winter, and thus any variation in the snowpack would result in relatively small changes of SCF, ultimately resulting in relatively small RSD. However, the snowpack that occurs during spring and autumn is subjected to stages of frequent melting and accumulation, and thus a thin snowpack, even with small levels of variation, could ultimately lead to relatively large RSD. From Fig. 6 we can also see that s2 shows good agreement with observations (Figs. 6a and c), apart from slightly larger RSD than observed over Central Canada (Fig. 6d). Meanwhile, the RSD given by g2 is less than the observation over all regions, revealing that g2 simulates much weaker interannual variability than s2.

4.2.4. Trends

Long-term changes in snow cover can be considered an indicator of climate change (Robinson and Frei, 2000). We assessed the models' performance in this regard (Brown, 2000; Dye, 2002), but concentrated on their abilities to simulate the variation trends of snow cover over the NH during winter months.

Figure 7c shows clearly that there is a significant decreasing trend in SCF over the region 30^°-60^°N during winter according to the monthly gridded NOAA SCF data for the period 1971-94. The most significant areas of decrease are located in Western Europe, the Tibetan Plateau, and eastern Northern America. Compared to the observation, s2 is able to simulate the decreasing trend of SCF well in these areas of decline (Fig. 7b), while g2 shows an unobvious decreasing trend over these areas (Fig. 7a). In addition, there is an obvious increasing trend of SCF in western North America-a trend that is barely simulated by s2. Although g2 is able to simulate the increasing trend in western North America, it is too weak compared to observations. In short, whether a decreasing or increasing trend is presented, g2 will show a weaker variation of that trend at the decadal scale compared to s2.

In general, both g2 and s2 are able to simulate the spatial pattern of SD and SCF well over the NH. The seasonal cycle of SD/SCF modeled by g2 is always higher than that by s2 over all months except the summer months (June-July-August; JJA). The root-mean-square error RMSE results reveal that the error of the seasonal cycle of SCF simulated by g2 is smaller than s2, and s2 has a smaller error than g2 in terms of the SD seasonal cycle. Finally, g2 produces weaker interannual variation compared to s2.

It is difficult to locate the sources of problems in GCMs. Even a very small modification of one parameter in a model can bring an undesired effect on another. Moreover, many parameters are needed in GCMs, and it is hard to acquire enough relevant observational data for their validation; thus, it is difficult to know which parameters are behaving improperly. In terms of the present study, discrepancies in SD/SCF between the two models and observations have been identified and described. In this section, we discuss the reasons that might lead to the differences between the two models.

Compared to s2, g2 reproduces excessive SD/SCF over all months except JJA in simulations of the seasonal cycle of SD/SCF. There may be several reasons for this, such as unreasonable simulation of precipitation, improper simulation of snow melting processes, and so on. However, g2 and s2 have the same snow processes but different atmosphere components. Therefore, the causes of the discrepancy between the two models should be explained by looking into their atmospheric forcings.

There is a significant difference between the two models' simulated SDs in winter from 1956 to 2005 over mid-high latitudes in the NH (Fig. 8a); the SD obtained by g2 is always higher than that given by s2. As a correct representation of snow cover requires accurate simulation of both precipitation and surface temperature, the variations in precipitation and temperature of the corresponding period over the mid-high latitudes of the NH were analyzed. Precipitation simulated by the two models shows an increasing trend, and the rate of increase according to g2 and s2 is 0.01 mm (10 yr)^-1 and 0.02 mm (10 yr)^-1, respectively. Although precipitation simulated by s2 is less than g2 in the early years, the difference by the end of the period is relatively small because s2 has the larger growth rate (Fig. 8b). For temperature, the values simulated by s2 are positive during the test period, while g2 produces negative values. The difference in temperature in the early years is about 2^°C, and then finally reaches about 3.5^°C because the rate of increase is 0.205^°C (10 yr)^-1 for g2 and 0.414^°C (10 yr)^-1 for s2 (Fig. 8c). When two models show little difference in precipitation, a lower temperature is good for converting more liquid water to snowfall and hinders the melting of ground snow. Consequently, more solid snow would accumulate on the ground, while a higher temperature would reduce the fraction of precipitation that falls as snow and speed up the melting of snow, and eventually lead to less snow accumulating on the ground. Thus, the temperature difference between the two models was considered as the main driving force behind a greater SD being simulated by g2 compared to s2. Furthermore, in the two models, SCF is completely dependent on SD, being given by where z _lnd is the surface roughness length and is set to a constant of 0.01 m in the two models. Obviously, in g2 and s2, the relationship between SCF and SD is monotonically increasing and highly nonlinear; SCF increases with increasing SD. When SD is relatively shallow, SCF increases quickly with increasing SD; however, when SD is relatively thick, an increase in SD can have a minor impact only on SCF. In the present study, SCF increases slowly with the increased SD; therefore, compared to s2, g2 also simulates much more SCF than SD.

The major physical difference between the two models is their cloud-related processes. A two-moment cloud microphysics scheme adopting different physical-based aerosol activation parameterizations is employed to estimate the aerosol indirect effects in GAMIL2, while in SAMIL2 the modified Edwards-Slingo scheme is used, which only considers the direct effect of aerosols. Aerosol indirect effects influence the radiation and surface energy balance by modifying radiative properties, as well as the amount and lifetime of cloud by changing cloud microphysical characteristics (IPCC, 2007). The difference between direct and indirect effects of aerosols could be reflected by the simulation of radiation flux. Figure 10 shows the difference in net radiation flux between the two models (s2 minus g2). The average net radiation difference over the regions 60^°-90^°N, 40^°-60^°N and 0^°-40^°N are 2.41 W m^-2, 0.25 W m^-2 and -3.65 W m^-2, respectively. In other words, compared to s2, g2 has a weaker heating effect on the land surface, which results in a lower surface temperature in mid-high latitudes of the NH.

In short, compared to s2, g2 has a weaker net radiation flux at the land surface in the mid-high latitudes of NH because it considers the indirect effects of aerosols, which results in lower simulated temperatures in g2. However, the precipitation difference between the two models is very small. Lower temperature is beneficial for converting more liquid water into snowfall, but hinders the melting of ground snow, and ultimately results in greater SD/SCF.

As reported in section 5.1, g2 reproduces more SD/SCF over all months except JJA compared to s2. However, compared to observations, g2 produces a smaller error than s2 in its simulation of the seasonal cycle of SCF, while s2 produces a smaller error than g2 in its simulation of the seasonal cycle of SD. In this section, we provide some explanations for these results.

The SCF is defined as the fraction of a model grid cell covered by snow, and is considered to be one of the largest sources of uncertainty in snow modeling (Niu and Yang, 2007). The variation of SCF can result in changes to surface albedo and significantly impacts upon the energy balance between atmosphere and land. In models, the land surface air temperature is calculated according to the surface energy balance equation. Therefore, whether or not the surface energy flux is well simulated could be reflected by the simulation of land air temperature. In the present study, through an analysis of global land temperature anomalies during the period 1956-2005, relative to the mean temperature of the years 1961-90 (Fig. 9), it can be seen that g2 is able to capture the variation of land air temperature well during the period 1956-2005, and the correlation coefficient between observation and simulation reaches 0.89. Meanwhile, s2 simulates a higher temperature compared to the observation, and the correlation coefficient is 0.85, which is slightly lower than g2. These results indicate a better performance of g2 compared to s2 in terms of simulating the surface energy balance. Moreover, as a small variation in SCF can impose a large impact on the surface energy balance, a good performance in simulating the annual cycle of SCF is achieved by g2.

Furthermore, due to the limited availability of SD data on the continental scale for validating the relationship between SCF and SD, great discrepancies exist among various SCF schemes. At present, most SCF schemes are defined as a function of SD and the ground roughness length, while a few have taken the variation of snow density and topography into account (Liston, 2004; Niu and Yang, 2007). (Li et al., 2009) showed that SCF obtained from Eq. (2) is less than observations, as well as that simulated by those schemes that consider the impact of snow density and topography, under the same SD. Therefore, the SCF simulated by g2 and s2 is less than observed for almost all mid- and high-latitude areas of the NH (Fig. 1). To acquire a more realistic simulation of SCF, the current SCF scheme in the two models should be further improved.

In most areas of the mid-high latitudes in the NH, both g2 and s2 produce more precipitation than observed (Table 1), and compared with s2, g2 has the smaller values (Fig. 8c) and better performance (Fig.9) in simulating land air temperature. Thus, much more precipitation will be converted to snowfall during winter/spring in g2 simulations, leading to more snow accumulating on the ground and higher SDs than observed. If we suppose that s2 simulates the same temperature as g2, it would also produce larger SD results. However, s2 presents significantly higher temperatures (Fig. 8c), which causes less precipitation to be converted into snowfall, and ultimately shallower SDs on the land surface. It is the offset of more precipitation and higher surface temperatures that lead to the good performance of s2 in simulating the annual cycle of SD, but it does not mean that s2 has the more accurate physical processes of SD accumulation.

Figure 9.  Global land temperature anomaly time series for g2 and s2 compared against the CRUTEM3 observational estimate. Anomalies are relative to the 1961-90 average for both models and observations. Shading shows the range across the ensemble. Units: ^°C.

Figure 10.  The difference in net radiation flux (units: W m^-2) between s2 and g2.

During 1956-2005, the global concentration of CO₂ increased largely from 312 ppm to 379 ppm. The response of the two models (g2 and s2) to increased CO₂ concentration is called "climate sensitivity", and is determined by the internal feedback processes that amplify or dampen the influence of radiative forcing on climate (IPCC, 2007).

(Chen et al., 2013) investigated the responses of global mean surface air temperature to idealized CO₂ forcing by using the output of abruptly quadrupling CO₂ experiments, and estimated the radiative forcing of quadrupled CO₂ and equilibrium sensitivity by adopting the Gregory-style regression method to understand the different climate sensitivities of g2 and s2. They found that the climate sensitivities of g2 and s2 are about 3.7 K and 4.5 K, respectively. It was clear that g2 has a weak response to this forcing and gains a little climate sensitivity. Through careful analysis, the authors concluded that it is the shortwave cloud feedback that contributes most of the uncertainties in the estimated net feedback. ENSO is also associated with some of the most pronounced year-to-year variability of climate features in many parts of the world (IPCC, 2007). According to (Li et al., 2013b), compared with the previous version of FGOALS (g1), g2 produces a weaker ENSO, e.g., the standard deviations of Nino3 index are 2.1 K and 0.83 K for g1 and g2, respectively. The weak response of g2 to the forcing makes the whole climate system change little from year to year. Therefore, g2 presents a weaker interannual variability of SCF. The inference is that the improper parameterization of cloud-related processes of g2 lead to the too-weak simulation of shortwave cloud feedback and the weak climate sensitivities, meaning more sensitivity tests are needed to further analyze and improve the model in this regard.

Based on the historical simulations of CMIP5, the SD and SCF simulated by g2 and s2 were validated using USAF/ETAC SD data, and SCF data from NOAA visible and MODIS products. The main conclusions can be summarized as follows:

(2) The spatial patterns of SD and SCF over most areas of the NH are simulated well by both g2 and s2, and the average spatial correlation coefficient over all months reaches about 0.7 and 0.8 for SD and SCF, respectively. However, the complex terrain and paucity of snow measurements on the Tibetan Plateau cause a large discrepancy between the two models and observations over this region.

(2) Compared to s2, g2 simulates a larger annual cycle of SD/SCF over the non-summer seasons, and the difference is most pronounced in late winter and spring. Analysis of RMSEs revealed that the seasonal cycle of SCF simulated by g2 contains less error than s2, while s2 produces a smaller error in its seasonal cycle simulation of SD. As g2 considers the indirect effects of aerosols in its atmospheric component, it produces a weaker net radiation flux at the land surface compared to s2, which results in lower temperatures simulated by g2. Meanwhile, the precipitation simulated by the two models are much the same, and thus the lower/higher temperature causes more/less liquid water to be converted into snowfall and hinders/accelerates the melting of ground snow, consequently causing more/less snow to accumulate on the ground.

(3) The analysis of interannual variability and long-term trend of SCF showed that both models simulate the significant decreasing trend well over (30^°-70^°N) in winter during the period 1971-94. However, as g2 has a weak response to an increasing of the CO₂ concentration and lower climate sensitivity, it presents a weaker interannual variation than s2.

The main purpose of the present study was to evaluate the fundamental variables of SD and SCF, which can be used to characterize snow cover. The role of atmospheric forcing on snow simulation was the main emphasis because the same land surface schemes are used in the two versions of FGOALS studied. In fact, an inaccurate description of snow cover itself has a great impact on temperature, the hydrological cycle and atmospheric circulations by affecting the surface energy and water balance. Improvements in snow parameterizations, including the parameterization of snow albedo, snow cover fraction, snow density, snow heat conductivity, the number of snow layers etc. are also necessary to improve snow simulations. In addition, although the snow-albedo feedback mechanism in coupled climate models is also very important, as the process is very complex, further research is required. Compared to the previous version, great improvements have been made for g2 in cloud-related processes and the model's climate sensitivity has been reduced. However, a too-weak response to CO₂ forcing still exists at present. Therefore, it is necessary to further improve the parameterization scheme of cloud-related processes to increase the response of the model to external forcing.