Study site

Sagarmatha national park is situated in the Solu–Khumbu district, in the northeastern region of Nepal (Fig. 1), and occupies the northernmost part of the Dudh Koshi river basin, which in turn is part of the Koshi river basin (or Sapta Koshi basin), one of the seven major hydrographic basins into which Nepal is subdivided. SNP, officially instituted on 19 July 1976, covers an area of 1141 km² and has unique geographical features, being surrounded on all sides by the highest mountain ranges on Earth. The terrain is extremely irregular, with altitudes ranging from 2845 m (Monjo) to 8848 m (Sagarmatha is the Nepalese name for Qomolangma (Mount Everest)).

Fig. 1. Location of the historical meteorological stations near Sagarmatha national park (SNP).

Reference ByersByers (2005) describes the climatic features of SNP: geographically, it lies within the subtropical Asian monsoon zone characterized by pronounced summer rainfall maxima, with over 80% of the annual precipitation falling during an approximately 4 month period between June and September (Mani, 1981; Reference Barry and ChorleyBarry and Chorley, 1982). Winters are normally dry, although occasional mid-latitude cyclones, driven by the subtropical jet stream that takes its winter position just south of the Himalayan ridge, can cause heavy snowfall events (Reference Klaus and HellmichKlaus 1966; Reference Zimmermann, Bischel and KeinholzZimmermann and others, 1986). A gradual transition from the dry winter season to the summer monsoon results from pre-monsoon convective precipitation events which are often accompanied by thunderstorms (Reference Ueno, Shiraiwa and YamadaUeno and others, 1993). The prevailing axis of the monsoons is southern and southwestern (Reference MüllerMüller, 1980). Reference AgetaAgeta (1976) observes that while monsoon activity from the southern foot of the Himalaya decreases, approaching the interior of the main range of the Himalaya precipitation increases at some of the higher points. For instance, the average annual precipitation of Namche Bazar (27.83° N, 86.72° E; 3450 m) is ∼1000 mm a⁻¹ and decreases with elevation (to ∼500 mm a⁻¹ in Dingboche (27.89° N, 86.83° E; 4355 m), while the total amount of precipitation around peaks and ridges can be four or five times greater than that in valleys (Reference Higuchi, Ageta, Yasunari and InoueHiguchi and others, 1982). Precipitation occurs mainly in the daytime on the ridges, but at night in the valleys.

The glaciers of SNP, like those of the Himalaya in general, exhibit certain morphological features of particular interest. Nearly all (28 out of a total of 29 in SNP) are ‘black glaciers’, known also as D-type or debris-covered glaciers; these are glaciers in which the ablation zone is almost entirely covered by surface debris that significantly alters the energy exchanges between the ice and the atmosphere. According to Reference MoribayashiMoribayashi (1974), Reference InoueInoue (1977) and Reference Smiraglia, Baudo, Tartari and MunawarSmiraglia (1998), these Himalayan glaciers are further characterized by a marked difference between their maximum elevation and the snowline elevation, a stable correlation between the accumulation area and the overall length of the glacier, abundant accumulation driven by wind and avalanches and multiple confluences that sometimes make it difficult to identify the principal lobe. Another characteristic feature is the debris cover which gradually increases toward the lower part of the glacier’s ablation zone. Here, the coverage becomes complete, and the glacier is regarded as inactive, presenting an apparently chaotic arrangement of detritus, depressions and cavities that often contain very small lakes of muddy water up to about 10 m in depth (Reference Iwata, Watanabe and FushimiIwata and others, 1980). Reference Fujii and HiguchiFujii and Higuchi (1977) report that the C-type glaciers (i.e. glaciers without debris cover on their ablation zone) of SNP at the end of the 1950s occupied an average area of 0.43 km² (standard deviation 0.55 km²). Because the minimum threshold for considering glacier areas in this study was set at 1 km², as we shall discuss below, the only debris-free glacier we consider is Langmuche glacier.

Cartographic supports and georeferencing

The topographic and thematic maps available for Nepal are well documented by T.B. Pradhananga (http://www.gisdevelopment.net/aars/acrs/2002/mgt), who provides a complete report of all the existing (topographical, land resource, thematic, derived, etc.) maps by Nepal’s Topographical Survey Branch. Of the existing cartography, only the Khumbu Himal map, published in 1978 and the current official map of Nepal published in 1997 are sufficiently representative of the area of SNP. Table 1 details the characteristics of these two maps, while Table 2 gives the nomenclature of the glaciers.

Both maps have the same scale, so it is possible to compare the level of detail of the information they contain. Compared to the Khumbu Himal map, the official map of Nepal provides slightly less definition of qualitative aspects such as the morphological representation of rocks and, of more importance for this study, no information on debris cover. This latter shortcoming means we cannot directly evaluate the variations over time in the areas of the accumulation zones, generally without debris cover, and of the ablation zones, where debris is present.

The Khumbu Himal map covers 98.0% of the park’s territory, extending from longitude 86°33′ to 87°07′ E, while the boundaries of SNP extend in longitude from 86°31′ to 87°01′ E. However, this map does cover the entire latitude range of the park. To take this into account and permit comparisons with the current official map of Nepal, the data referred to the latter map have also been restricted, in the analysis and the tables that follow, to the longitudinal limit (86°33′ E) of the Khumbu Himal map.

In order to compare the glacier areas it is first necessary to define the time periods to which the maps refer. As Table 1 shows, the Khumbu Himal map takes in a wide span of time between the photographic survey (1921), the terrestrial photogrammetric survey (1935 and 1939) and the terrestrial survey, itself done over a lengthy period, 1955–63 (Reference Schneider and HellmichSchneider, 1967). To narrow down the dating of this map, we use the fact that Müller, as we describe below, attributes the data of his inventory conducted on the Khumbu Himal map (Reference Schneider and HellmichSchneider, 1967) and the subsequent analyses to ‘the 1950s and early 1960s’ (Reference MüllerMüller,1980). On the strength of the period given by Müller, we consider the Khumbu Himal map to be representative of the late 1950s (and hereinafter refer to it as the 1950s map). The period of the current official map of Nepal is better defined, and can be considered representative of the beginning of the 1990s (and is hereinafter referred to as the 1990s map). A comparison of the two maps thus allows us to analyze the changes occurring over a period of nearly 35 years.

The two maps were georeferenced according to the projection system of the official map of Nepal. However this was not done by simply altering the projection system of the Khumbu Himal map. In fact, J. Gspurning and others (http://avalanchemapping.org/Repapers/GISinhighmountainKumbu.pdf) have pointed out the horizontal deviations which exist between identical map control points, due to the poor accuracy of the local projection systems. We therefore proceeded to resample the Khumbu Himal map using the ArcView^© tool, based on >100 map control points (the mountain peaks) evenly distributed on the maps. The a posteriori check showed that the maximum horizontal deviation (along the x and y axes) did not exceed 15 m (georeferencing error (GE)), compared with 336 m for the y axis and 19 m for the x axis observed by Gspurning and others.

Evaluation of the accuracy of the maps

Before analyzing the precision errors arising from the process of cartographic representation and its subsequent conversion from analogue to digital form, we must consider the problem of the interpretation of glacial features (accuracy). This factor must be taken into account to ensure that any differences observed between the two maps effectively reflect real changes, rather than artefacts of different degrees of accuracy in the construction of the two maps. For the specific purposes of this work, this implies evaluating the comparability of the glacier outlines on the two maps. This is unfortunately not a simple matter, as it would require retracing the process by which the cartographers discriminated between what is and what is not mapped as a glacial feature. Firstly, there are a multiplicity of data sources to consider: photographic material, photogrammetric images, direct observations, as well as the personal experience of the cartographer in interpreting the features on the terrain. Secondly, any assessment would have to include various forms of ice cover in order to be representative of the entire park area. Finally, we have to consider that the photogram-metric images used for constructing the Khumbu Himal map were taken on glass plates whose current state of preservation does not allow them to be easily re-examined.

In metrology, the interpretation of glacial features on maps is commonly evaluated in terms of a data-accuracy index, defined as the discrepancy between the measured value (mapped glacier outlines) and the real value of the studied quantity (actual glacier outlines). Such a discrepancy generally arises from biases introduced by the instrumental method or, as in our case, by human factors (interpretation of glacier features). The error is considered systematic because it remains constant on repeating the measurement, making it difficult to determine not only its magnitude but even its very existence, except through a repetition of the entire measurement process (which is not feasible in our case) or a consistency evaluation of the results (Reference Glosup and AxelrodGlosup and Axelrod, 1998). Consequently, despite the difficulty of directly analyzing how the glacial features were interpreted on the two maps, it is important to determine whether the results of the comparison can be ascribed to real changes, so in the discussion of the results we attempt to evaluate the congruence between the results and the processes that may have brought them about.

Calculation of the cartographic error

A consideration of fundamental importance for Geographical Information System (GIS) processing is the precision of data. In the present work, the data processing involves conversion of topographic maps from analogue to digital form, entailing a series of approximations, with inevitable effects on data quality. A final cartographical representation is dependent on various factors, including the precision of the map, the georeferencing error and the precision of the operator, which gives rise to the vectorization error. Assuming the vectorization error to be minimal, and having dealt with the GE of the maps as discussed above (15 m along the x and y axes), we can now further account for the precision of the data.

Every map has a well-defined metric precision which derives from the process by which it was constructed. The precision of the location of points and the level of detail of a topographic map, whether in digital or analogue format, are determined by its scale factor. The precision of a derived cartographic representation is the degree of planimetric and altimetric tolerance between the coordinates of any given point on the source map and the coordinates of that same point on the derived map. Evaluating the precision of a derived map is therefore necessarily dependent on an evaluation of the precision of the source map. The graphical resolution limit of a map is assigned an arbitrary value of ∼0.2 mm, based on the threshold visible to the human eye. The level of approximation of a map is conventionally given by this limit of 0.2 mm multiplied by the scale factor. Therefore, for the source maps used in this work (scale 1 : 50 000), the approximation, i.e. the linear resolution error (LRE), corresponds to 0.2 mm × 50 000 = 10 m (Reference InghilleriInghilleri, 1974).

The overall linear error (LE) associated with the distance between two points must take into account the GE discussed earlier, as well as the LRE, and is the root mean square (rms) of the two errors (Reference Bevington and RobinsonBevington and Robinson, 1992):

GIS-based processing is often used for analyzing variations in areas and volumes. It is necessary to determine the deviation between the actual values and those computed based on the digitized data. Because the perimeter of an area is delineated by manual digitizing, the value of the planimetric vectorization error is often used to quantify the areal resolution error (ARE), which indicates the range of tolerance of the area measurements in question. The ARE is usually computed as the product of the perimeter, l, of the area in question and the LRE (Reference InghilleriInghilleri, 1974). Since the LE is greater than the LRE, the areal error (AE) associated with an area must therefore be computed as:

In Table 3, the AE is referred to the glacier area variations between the end of the 1950s and the early 1990s, computed as the rms of the AEs associated with the 1950s and 1990s areas. The overall error was also referred to the 1950s area of each glacier, thus directly compared with the glacier area difference.

The digital elevation model

The digital elevation model (DEM) was constructed by digitizing the contours of the 1990s map. The interpolation was done by the kriging method using the Surfer^© software, setting the north–south isotropy of the data identified from analysis of the semivariogram. A DEM pixel dimension of 20 m × 20 m was used. Considering that the mean elevation range of the glaciers is ∼1620 m for a mean length of ∼6.7 km, there were on average ∼8 cells interpolated between adjacent contour lines. To judge the validity of the above processing method and hence the precision of the morphometric data thus obtained, we computed a posteriori the average height deviation (along the z axis) between the interpolated values and ∼100 map control points, obtaining an average absolute difference of 6.3 m. The mean error associated with the average glacier elevation (Table 2) is therefore 0.2%. For the surface gradient, because the calculation is based on the relative values between consecutive DEM pixels, the accuracy of the DEM absolute values is less important (Reference QuinceyQuincey and others, 2007) and has here been disregarded.

Statistical analysis

The statistical method used in this work is a correlation analysis whose results are expressed as a coefficient of determination, R ². The frequency distributions of the dependent variables are also given to show that all correlations are based on Gaussian distributions. The correlations found are validated by a parametric inferential analysis of the observed correlation coefficient, using a Student’s t test (Reference Venables and RipleyVenables and Ripley, 2002). This test evaluates the probability of the null hypothesis H₀: r = 0 (i.e. absence of linear correlation between the two examined variables, meaning that the correlations found are due to chance), and hence the probability that the correlations found are significant (H₁: r ≠ 0). A significance level (α) is generally chosen to define the rejection range of the null hypothesis as α = 0.05, α = 0. 01 or α = 0.001. The p value gives a measure of the plausibility of the null hypothesis, defined as the minimum significance level, α, of the test for which the null hypothesis would be rejected. If p < α the null hypothesis is rejected and the correlation is deemed significant. The significance of the correlation coefficient, r, is obtained using the table of critical values of r_{α ,n} , as a function of the number of observations, n, and the different significance levels, α. For the linear correlations discussed below, the observed correlation coefficient, r, was compared with its critical values, r_{α ,n.} If r> r_{α ,n} then p < α, and the null hypothesis was rejected. All the graphs included in this work were found to be statistically significant, and in each case we give the value of the coefficient of determination, R ², and the significance level of the associated correlation coefficient (p < α).

The meteorological data

There are not many weather stations above 2000 m in Nepal, as noted by Reference Shrestha, Wake, Mayewski and DibbShrestha and others (1999). They provide a list of stations operated by the Department of Hydrology and Metereology (DHM), some of which have been active since the mid-1960s. Temperature and precipitation data are available as monthly means (Nepal: DHM, 1998). Figure 1 shows the locations of stations managed by the DHM situated near SNP. Also included in the figure is Shegar station (28.63° N, 87.08° E; 4300 m a.s.l.) in Tibet, for which annual temperature and precipitation graphs are described by Reference Jin, Xin, Che, Wu and MoolJin and others (2005). In the discussion we present annual mean temperature and precipitation graphs for Kathmandu station, as this recorder is the meteorological reference for Nepal. In this case, after a 9 year period of joint operation, a station located at the Indian Embassy (IE) (27.44° N, 85.20°E; 1324 m a.s.l.) was replaced by one located at Kathmandu International Airport (IA) (27.42° N, 85.22° E; 1336 m a.s.l.). We developed a linear regression model of the annual mean temperature and annual total precipitation data for the IE station, based on the IA annual data for the overlap period (1968–75). The extension of the IE record using the regression model (Kathmandu-simulated series) allows us to estimate the annual temperature and precipitation trends during the observation period. The error arising from the estimates is the rms of the differences between the reconstructed values and the observed values during the overlap period (Reference Hays and WinklerHays and Winkler, 1970).

Glacier identification

Previous work has constructed inventories of the glaciers in SNP. Reference MüllerMüller (1970) reports an inventory of glaciers in the Qomolangma region (corresponding in extent to the current SNP) which can be considered one of the pilot studies in this field. This work was undertaken using the map of the late 1950s on a 1 : 50 000 scale (Reference Schneider and HellmichSchneider, 1967). Reference Fujii and HiguchiFujii and Higuchi (1977), referring to this inventory, analyzed the relations between glacier formations and the topographic conditions of the glacier basins. Unfortunately for the purposes of the present study, it was not possible to use the area data from the above inventory directly, because it restricts the area calculation to those glaciers classified as ‘extensive’, giving only the geographical position for the smaller glaciers. Such a method results in a substantial underestimation of the total SNP glacier area (287 km²). In 2001 the International Centre for Integrated Mountain Development (ICIMOD; Reference Mool, Bajracharya and JoshiMool and others, 2001) constructed a glacier inventory for all Nepal, and thus also for SNP. Working at the national level, a variety of data sources were used including satellite images, aerial photographs and maps. The total SNP glacier area was calculated to be 360 km². However, this inventory also proved unsuitable as a basis of area comparison for the present work, because the multiplicity of its sources cannot be attributed to a precise historical period. That said, since it is not the purpose of this work to construct a new inventory, but only to evaluate changes in glacier areas, it was decided to adopt the glacier codes proposed by ICIMOD and shown in Table 2. We see in Table 2 that glacier nomenclature differs considerably across the two maps. We choose to use the ICIMOD inventory as a reference, giving it preference in assigning unambiguous place names. For those glaciers that were not assigned a name in the ICIMOD inventory, we used, if reported, the name given in the maps (e.g. Phunki, Machhermo, Kyajo, Chhuitingpo) or the ICIMOD glacier code (e.g. Kdu_gr 181, Kdu_gr 125, Kdu_gr 38). The only remaining problem with the ICIMOD inventory concerned Ama Dablam glacier, which was swapped with Ombigaichain glacier, and has here been called Duwo glacier in accordance with both maps. Table 2 also shows the coordinates of the glacier centroids as identified on the 1990s map. Figure 2 shows the glacier areas at the end of the 1950s and in the early 1990s.

Fig. 2. Glacier comparison map of SNP.

For classification criteria used to define the extent of each glacier, we referred to the glacier classification guidance for the GLIMS glacier inventory (Reference Rau, Mauz, Vogt, Khalsa and RaupRau and others, 2005). This gives different classification methods, depending on the purpose. For analyzing the extent of glaciers, the ‘Form parameter group’, which essentially describes the outline of the glacier, is suggested. Following the criteria given for this category, we associated glaciers with their catchment areas, thereby defining compound basins (fig. 20 of Reference Rau, Mauz, Vogt, Khalsa and RaupRau and others, 2005) in which several individual accumulation basins feed into a single glacier system. In accordance with Reference Mayer, Lambrecht, Beló, Smiraglia and DiolaiutiMayer and others (2006), the tributary valley glaciers are included as part of the main glacier even if they are not directly connected to the accumulation and ablation zone of the main body in one of the two maps. Proceeding in this way, these areas were considered to be the glacier’s possible ‘feeding’ zone. Using the above classification criteria also ensures that the comparison is being made between the same glacier surfaces, without mistaking the joining or detaching of tributary glaciers as increases or decreases in glacier area. For instance, on the basis of these criteria, Gyubanare glacier, which is usually considered independent (Reference MüllerMüller, 1970; Reference Fujii and HiguchiFujii and Higuchi, 1977), is not detached from Ngojumba glacier.

Working on the early-1990s map, we began by digitizing the surface hydrography of SNP and subsequently the corresponding hydrographic basins. We then proceeded to digitize the glacier areas on both the 1950s and 1990s maps. As discussed above, to permit a comparison between individual glaciers on the two maps, it was necessary to identify glacier basins, choosing as the basin outlet the glacier front which had advanced the most between the two periods under study. By so doing it was possible to ensure that each basin entirely contained the area of the corresponding glacier in both maps. To make comparisons between individual glaciers, we chose to consider only those glaciers whose area in the late-1950s map exceeded a threshold of 1 km². Proceeding in this way, a total of 29 glaciers, accounting for nearly the totality of glacier area in SNP (98%), were analyzed.

In Table 3 we see that the largest glacier was Ngojumba glacier, with a 1950s area of 98.7 km². In the 1950s the average elevation of the glacier areas was 5488 m, with Chhule glacier having the lowest (5099 m) and Khumbu glacier the highest average elevation (6100 m). The average aspect of the glaciers was 172°; in other words the glaciers were, on average, oriented southwards. Specifically, Tingbo and Cholotse glaciers were exposed to the west, Cholo glacier to the east and Kdu_gr 38 glacier to the northeast, while the remaining majority of glacier surfaces were instead exposed to the southwest. The average glacier slope during the 1950s was 28.0%. The two extremes were Kyajo glacier, with the lowest average slope (18.8%), and Phunki glacier, with the highest average slope (41.7%).

Table 3 also shows the morphometric characteristics of the glaciers in the 1990s. We observe that the total glacier area in the 1990s (without the 1 km² cut-off threshold) was 384.6 km² (AE = 20.7 km²), and the total glacier area beneath the threshold was 15.7 km²; the average glacier (above the threshold) area was 13.5 km², but with standard deviation (20.0 km²) nearly twice the mean. There is considerable variability in the size of individual glaciers, which can differ by nearly two orders of magnitude. In the 1990s, Ngojumba glacier was still the largest, with an area of 104.7 km². The average glacier surface elevation was 5491 m, with Cholo glacier having the lowest (5002 m) and Khumbu glacier the highest average elevation (6154 m). The average aspect computed in the 1990s was 176°. There were eight glaciers exposed to the south, six to the west, five to the southeast and five to the southwest. Langmuche, Langdak and Kdu_gr 38 glaciers were exposed to the northeast and Cholo and Thyangbo to the east. The average slope of the glaciers in the 1990s was 24.8%, with the two extremes being Kyajo glacier, with the lowest average slope (10.5%), and Phunki glacier, with the highest average slope (41.2%).

Glacier area variations

Figure 2 shows a visual comparison of glacier areas, while Table 3 gives the percentage change in the area of each glacier from the 1950s to the 1990s, and the morphometric parameters for these two periods. We observe that:

In summary, our comparisons revealed that during the study period SNP experienced a small net reduction in glacier cover. The changes to individual glacier areas were variable, with a reduction in size of a large number of glaciers being offset by an increase in size of a smaller number of glaciers, thus mitigating the overall loss of glacier area. In the following subsection, we examine these opposing variations in more depth, attempting to correlate them with the morphological characteristics obtained from the DEM.

For each glacier, the change in surface area from the 1990s to the 1950s is compared with its size (area) in the 1950s. Figure 3a shows the percentage change in glacier area as a function of glacier area. For the smaller glaciers (<20 km²) there is a clear tendency toward loss of glacier area over the study period, but no correlation between glacier size and area change is evident. The larger glaciers (>20 km²), however, are less likely to have lost area over the study period, and overall there is no clear trend toward loss or gain. In Figure 3b we see that the average area loss of the glaciers was 18.2%, although the effect of the larger-sized glaciers, which on average increased in area, meant that overall the park lost only 4.9% of its glacier area.

Fig. 3. (a) Glacier area variation from the 1950s to the 1990s, as a function of glacier size. (b) Frequency distribution of glacier area variation.

Figure 4a shows the percentage change in glacier area over the study period, as a function of glacier average elevation. Lower-altitude glaciers lost the most area, while higher-altitude glaciers increased their area (R ² = 0.34; p < 0.001). The glaciers of SNP had an average elevation of 5488 m in the 1950s, ranging from 5099 m (Chhule glacier) to 6100 m (Khumbu glacier) (Fig. 4b).

Fig. 4. (a) Glacier area variation from the 1950s to the 1990s, as a function of glacier elevation. (b) Frequency distribution of glacier mean elevations.

From Figure 5a we see an inverse relation (R ² = 0.21; p < 0.01) between percentage change in glacier area and average glacier slope. The glaciers which lost most area were those of highest slope, while the glaciers whose area increased were those of lower slope. The average slope of the glaciers was 28.6% in the 1950s, and the distribution is close to normal, but with three glaciers with slope >40% (Kdu_gr 38, Phunki and Tingbo glaciers) (Fig. 5b).

Fig. 5. (a) Glacier area variation from the 1950s to the 1990s, as a function of glacier slope. (b) Frequency distribution of glacier average slopes.

We also analyzed the relationship between glacier area and slope. From this comparison, it emerged that although there is, as expected, an inverse relation, its magnitude is not significant due to the presence in the park of several small glaciers with low slope gradients.

Figure 6a shows the glacier area in the 1950s relative to glacier aspect. The plot shows the larger glaciers were prevalently oriented to the south. The orientation of glaciers whose area decreased (empty circles) and increased (black squares) between the 1950s and the 1990s is indicated. Glaciers that increased are prevalently oriented toward the south. Figure 6b and Table 4 highlight this distribution, showing the area variation of each glacier (absolute value in km²). The glaciers which expanded were oriented to the south, while those which contracted were oriented toward other sectors. The greatest reductions were observed in the southwest sector, and amounted to 10.4 km², representing 39.0% of the total glacier area losses.

Fig. 6. (a) Glacier size as a function of glacier aspect. The empty circles represent glaciers that experienced area loss between the 1950s and the 1990s, while the black squares represent the glaciers whose surface area increased. (b) Grouping of total glacier area variations (km², radial axis) according to their orientation sector (°, circular axis).

Glacier area variations and changes in morphological features

Figure 7a shows the relationship between change in glacier area and change in glacier altitude during the study period. The glaciers which lost most area were, in the 1990s, situated at a lower average altitude than they were in the 1950s, while the glaciers which increased in area also exhibited an increase in average altitude (R ² = 0.29; p < 0.01). The area variations chiefly affected the accumulation zones of the glaciers rather than the ablation zones; in fact if the expanding glaciers increased their average altitude, they must necessarily have expanded more in their accumulation zones than in their ablation zones. Conversely, the contracting glaciers, whose average altitude decreased, must necessarily have suffered greater losses in their accumulation zones than in their ablation zones. We note, in any case, that the average elevation change (−1 m) is very small (Fig. 7b).

Fig. 7. (a) Glacier area variation from the 1950s to the 1990s, as a function of change in average elevation. (b) Frequency distribution of glacier elevation differences between the 1990s and 1950s.

Very similar behaviour is shown in Figure 8a, a plot of glacier area change over the study period against glacier slope change. Glaciers which shrank the most also showed a reduction in slope gradient, while the glaciers that expanded exhibited increased slope (R ² = 0.46; p < 0.001). The mean change was a reduction in slope of 0.034 (Fig. 8b).

Fig. 8. (a) Glacier area variation from the 1950s to the 1990s, as a function of change in slope. (b) Frequency distribution of glacier slope differences between the 1990s and 1950s.

These last two findings support the hypothesis that glacier area change, whether positive or negative, is primarily the result of changes occurring in the accumulation zones rather than in the ablation zones. This allows us to overcome the limitations of the cartographic representations. (As noted previously, the official map of Nepal, unlike the Khumbu Himal map, does not record any distinction between ablation zones (generally debris-covered) and accumulation zones (generally free of debris).)

The hypothesis that the greatest glacier area changes occur in the accumulation zones does not rule out the possibility that notable variations may also occur at the glacier front. Figure 9b is a map of Ngojumba glacier, the largest glacier in SNP, which expanded from 98.7 km² in the 1950s to 104.7 km² in the 1990s (+6.1 %, AE = 3.4%). Figure 2 shows its location in the north of SNP and its south-facing aspect (182°). We see that the area of the accumulation zone has increased. Figure 9a is a map of Cholo glacier, one of the smallest in SNP, which shrank from 2.1 km² in the 1950s to 1.5 km² in the 1990s (−31.2%, AE = 10.1 %). Figure 2 shows its location in the centre of the park and its east-facing aspect (83°). More specifically, Cholo glacier is situated in a western valley of the Khumbu valley and drains into the Imja Khola river. We observe here that the area of the accumulation zone has substantially decreased.

Fig. 9. Two examples of glaciers in which the greatest positive and negative area variations occurred primarily in the accumulation zones: (a) Cholo glacier, (b) Ngojumba glacier.

Divergence between individual glacier area variations

The glaciers of larger size, which were situated at higher altitudes, increased in size during the study period. Conversely, the smaller glaciers, situated at lower altitudes, showed a reduction in size (Figs 3 and 4). Similar findings by Reference Paul, Kääb, Maisch, Kellenberger and HaeberliPaul and others (2004) on a sample of 938 alpine glaciers from 1973 to 1998/99 and Reference Jin, Xin, Che, Wu and MoolJin and others (2005) for the Tibetan Plateau (∼900 glaciers) show a larger contribution by small glaciers to the total loss of area during the study period. The higher sensitivity of small glaciers to climate effects can be explained by considering that when the equilibrium-line altitude increases, it can extend above the altitude of many of the small glaciers, whereas the debris on the tongues of large glaciers often prevents them melting (Reference Jin, Xin, Che, Wu and MoolJin and others, 2005). This result indicates that the interpretation of the glacier outlines on the two maps can be considered realistic, as it is consistent with these other studies.

We observed that area increase was most apparent for the glaciers oriented toward the south (Fig. 6b) and occurred principally in the glacier accumulation zones (Figs 7 and 8). From Figures 2 and 5a we noted also that the larger glaciers (Ngojumba, +6.1%, AE 3.4%; Bothe Kosi, +12.6%, AE 3.4%; Lumsamba, +7.6%, AE 4.4%; and Nuptse, +25.4%, AE 7.3%), i.e. those exhibiting increased area, are predominantly oriented to the south. Conversely, area reductions during our study period were most evident in the glaciers oriented in other directions. Also in these cases, the area of the glaciers most affected was the accumulation zone. We need to determine whether the glacial features have been interpreted correctly on the two maps. Although there are a few studies (e.g. Reference HewittHewitt, 2005) which report variations in accumulation zones that exceed those in ablation zones, these can plausibly be accounted for in terms of climate change, as we show below. The fact that the glaciers whose area decreased showed the largest changes in their accumulation zones is consistent with their steeper slopes and lack of debris cover. Melting induced by temperature rises may have caused more extensive decoupling, fragmentation and avalanches here than in the flatter, debris-covered ablation zones, which are more susceptible to down-wasting than to area loss.

The fact that the larger glaciers are primarily oriented to the south (Fig. 6a) is accounted for by the southern orientation of the main valleys in the area, and the predominantly south to southwesterly origin of precipitation patterns (Reference MüllerMüller, 1980). The monsoon travels south–north in this area, before veering to the west (Reference BarryBarry, 1981). In addition, these glaciers are not occluded to the south by mountain ranges perpendicular to the monsoon that might act as a barrier against the mass of moist air which produces precipitation. If these glaciers have attained their present large size thanks to a favourable orientation with respect to monsoon precipitation, we conjecture that their increase in area may likewise have been favoured by an increase in precipitation. However, there was no increase in the area of glaciers situated at lower latitudes which, in addition to being smaller, have an orientation that diverges from the south–north axis. As shown in Figure 2, there are many glaciers which exhibit these characteristics (Chhule glacier, Nare glacier, etc.).

The above observations suggest that alongside temperature variations which are often considered to be the main cause of glacier change, variations in precipitation may also play a determining role. Temperature variations alone cannot account for why some glaciers increased in size while others (all of the same type) decreased in size. This latter observation is better explained by considering that there was also a precipitation increase over the study period. Figure 11 shows the annual precipitation trends for the stations situated near SNP and the Kathmandu-simulated series from the end of the 1950s to the early 1990s.

Fig. 11. Trend of total annual precipitation for the stations situated near SNP and for the Kathmandu series obtained as described in the text. For Shegar station, the graph was produced by Reference Jin, Xin, Che, Wu and MoolJin and others (2005). Data for the other stations are from Nepal: DHM (1998) and were processed for this work.

All the stations show an upward trend in precipitation, though for the series closest to SNP there is no complete dataset available. Table 7 shows, for each station, the yearly total precipitation for the periods for which data are available, and the corresponding observed increase per year over the period. It is not possible to make a direct comparison between these data, which come from different observation periods. We observe from Figure 11, however, that all the stations recorded a more pronounced increase until the mid-1970s, and a less steep rise in the following period. Here again, however, there remains uncertainty arising from the absence of data for the higher altitudes.

Glacier orientation relative to the prevailing direction of the monsoons may cause some glaciers to receive more precipitation than others, favouring those exposed to the south. What is more, increased precipitation may explain the increase in the area of the glacier accumulation zones. For these glaciers, the effects of the precipitation may predominate over effects of temperature, and the high elevations of these accumulation zones protect them from melting.

For the glaciers oriented in other directions and at lower elevations, which experienced a reduction in accumulation-zone area over our study period, the effects of temperature will have prevailed. Reference Fukui, Fujii, Ageta and AsahiFukui and others (2007) show, testifying to the effect of temperature at these altitudes, that in the Khumbu glacier valley the permafrost lower limit has risen 100–300 m between 1973 and 1991, stabilizing at 5400–5500 m. These glaciers, in addition, are often very steep, with notable discontinuities. As Reference MüllerMüller (1980) reported, this is a common feature of many glaciers in SNP: the extremely steep and large headwalls and side-walls of many valley glaciers contribute large amounts of avalanche snow to the accumulation area, frequently causing very irregular surface contours and, thus, irregular firn lines. An examination of Cholo glacier in Figure 9a shows a very steep glacier headwall in the 1990s, with remnant ice around the ridge-line at the highest elevation, and a disconnected lower-elevation glacier tongue, presumably fed by avalanches. For these glaciers, as explained above, it is therefore plausible that the accumulation zones, with steep slopes, may have lost more snow and glacial cover than the flatter, debris-covered glacier tongues, which instead have experienced greater loss in thickness (down-wasting) than in surface area (Reference Nakawo, Yabuki and SakaiNakawo and others, 1999).

A similar interpretation to that proposed above has been suggested in work relating to another period, which investigated the expansion of the glaciers of the Karakoram (Reference Diolaiuti, Pecci and SmiragliaDiolaiuti and others, 2003; Reference HewittHewitt, 2005). For these glaciers, situated at high altitude, an increase in winter precipitation was observed. Effects of the rise in winter temperatures resulting from global warming were found to be less critical than winter precipitation increase for these high-altitude glaciers (Reference Archer and FowlerArcher and Fowler, 2004). The glaciers at intermediate altitude, according to Reference HewittHewitt (2005), did not benefit from this expansion and instead continued their process of retreat.