Page 1 of 57
Projection of permafrost, its active surface layer and active layer
affected by climate change in the Swiss Alps.
J.E.D. Kerkmeijer (10753842)
Supervisor: Mr Dr J. H. van Boxel (MSc) Amsterdam, 27 June 2018
Page 2 of 57
Abstract
The disappearance of permafrost is an indicator of global climate change. Climate change is bringing along all kind of hazards, like mass movement originating from permafrost in the susceptible dense Alps. The aim of this paper is to provide a projection on how the distribution and depth of permafrost, the active surface layer and the active layer will be affected by climate change. The distribution, and depth of permafrost, active surface layer and active layer throughout the Swiss Alps was calculated starting from general equations that described heat transport in the soil. Of the 62 site points throughout the Swiss Alps, 42 had permafrost in 2018. If the mean temperature would rise 2 °C this would decline to 21 site points for the year 2100. In 2018, the occurrence of permafrost between the range of 2330 m and 2660 m elevation also depended on other variables than elevation only, like slope angle, aspect or surface temperature amplitude. This range shifted to 2536 m and 3100 m elevation for the year 2100. The active layer became thicker for every site over time, while the active surface layer did not change in depth. The influence of climate change on the active surface layer could not be well predicted, given that the temperature amplitude was assumed not to change over time. In 2100, permafrost will vanish below 2536 m elevation. In addition, where permafrost still occurs, the active layer will become thicker.
Page 3 of 57
Content table
Abstract ... 2 Content table ... 3 1 Introduction ... 4 2 Framework ... 5 2.1 Definitions ... 5 2.2 Nomenclature ... 62.3 Depth permafrost base ... 7
2.4 Homogenous soil ... 7
2.5 Depth active surface layer ... 8
2.6 Depth permafrost table ... 9
3 Methods ... 9
3.1 Descriptive ... 9
3.2 Predicting air temperature ... 9
3.3 Average surface temperature ... 10
3.4 Temperature amplitude ... 10
3.5 Calculations ... 10
4 Results ... 12
4.1 Descriptive ... 12
4.2 Predicting average surface temperature ... 12
4.3 Temperature amplitude ... 13 4.4 Calculations ... 14 5 Discussion ... 16 6 Conclusion ... 19 7 References ... 19 Appendix ... 22 1 extra graphs ... 22 2 spreadsheet of calculations ... 22 3 matlab script ... 22
Page 4 of 57
1 Introduction
Permafrost is defined as a ground layer that remains at or below 0 oC during two or more consecutive years (Osterkamp & Romanovsky, 1999). It lies under about a quarter of the land surface on Earth (Nelson, Anisimov & Shiklomanov, 2002) and since the last little ice age (during the mid-17th century) it has been degrading, caused by an increase in temperature (IPCC, 2013). Permafrost can be found in the tundra and in the alpine climate (Wubbe, 1987). Therefore is permafrost also present in the European Alps, whereby it is likely that Switzerland is the country in the Alps with the largest permafrost area followed by Italy, Austria, France, and Germany (Boeckli et al., 2012). Furthermore, Slovenia and Lichtenstein may have marginal permafrost areas (Boeckli et al., 2012). Additionally, like every soil, permafrost has an active surface layer, where the temperature fluctuates in the topsoil caused by daily and annual soil heat flux (Anisimov, et al., 1997). As soil heat flux is the energy exchange between the earth surface and the soil (Su, 2002).
Consequently, at high altitudes, the Alps’ air temperatures could rise above 0 oC especially during the summer (Gruber et al., 2004). Then, it is plausible that the permafrost has an active layer as well. As an active layer is a surface layer on the permafrost which means the temperature in that specific layer rises above 0 oC sometime during the year (Scherler et al., 2010), resulting in the freezing and thawing of the active layer (Osteramkp & Romanovsky, 1999).
The active layer and the permafrost thaw are considered to be important mechanisms that control the slope instability, landscape evolution and natural hazard potential in mountain areas (Gruber & Haeberli, 2007). As the Alps are densely populated for the most part (Gobiet et al., 2014), the area is susceptible to climate-related hazards (Gruber & Haeberli, 2007). In addition, the IPCC fifth assessment report summary for policymakers (2013) strongly suggests that permafrost's
temperature is increasing as for the thickness of the active layer, whereupon permafrost could vanish in many places as a result of global warming (Lunardini, 1996). This already resulted in an increasing number of rock slope failure events that the Alps have recorded over the past decades (IPCC, 2013). For example, during the remarkably dry and hot summer in 2003 in the European Alps, many rockfall events that originated from permafrost in steep bedrock, occurred without heavy precipitation or earthquakes (Gruber et al., 2004; Gruber & Haeberli, 2007).
Moreover, with earthquakes or heavy rainfall events, when water cannot infiltrate the soil, the magnitude of impact could be bigger by resulting mass movements like landslides and mudflows (Kääb, 2008; Crozier, 1999).
The European Alps have a long tradition of scientific research on landscapes by perennial surface ice and glaciers (Haeberli, 2005) as well as its influence on various natural and socio-economic sectors (Gobiet et al., 2014). In the European Alps, the amount of monitoring of the
permafrost has increased through the last decades due to its close relationship with climate change and the increasing interest in climate change (Haeberli, 2005). Still, there is an uncertainty in distribution models and a lot of uncertainty in future projections of the Alps region (Gobiet et al., 2013). Most models are founded on empirical statistical relationships, and provide an indication of permafrost with restricted accuracy and lack of knowledge about the depth of the permafrost layers (Boeckli et al., 2012).
The aim of this paper is to provide a projection of how the distribution and depth of
permafrost, the active layer, and the active surface layer will be affected by climate change. First, this paper provides some definitions that are assumed in this paper. Then, a representative temperature gradient through the soil is made. Afterward, formulas are created to give an output that represents the thickness of the active surface layer, the thickness of the permafrost table and the depth of the
Page 5 of 57
average surface temperature is increased. Finally, the simplification and assumption in the calculations are discussed before providing a conclusion.
2 Framework
2.1 Definitions
The definitions that are displayed in Table 2 were used to delimit this paper and to provide an overview for the reader.
Term Definition
Permafrost Permafrost is defined as a ground layer that remains at or below 0 oC during two or more consecutive years (Osterkamp & Romanovsky, 1999). Between permafrost table and base of permafrost in Figure 1. Temperature of 2100 The temperature will increase with 2.0 oC according to the climate
medium scenario (Gobiet et al., 2014).
Active surface layer Top layer where temperature fluctuates in the soil caused by daily and annual soil heat flux (Anisimo et al., 1997). Between the surface and upper limit of seasonally invariant temperature in Figure 1.
Soil heat flux is the energy exchange between earth surface and soil (Su, 2002). Active layer A layer between the surface and the permafrost where the temperature
rises above 0 oC at some time during the year (Scherler et al., 2010), also seen in Figure 1.
Permafrost table As Figure 1 shows it is the upper boundary of permafrost and the lower boundary of the active layer. Therefore, it is the dividing line between permafrost and possible active layer above the permafrost.
Page 6 of 57 2.2 Nomenclature
An overview of the symbols and their meaning are provided.
T = temperature [oC]
z = depth [m]
x = thermal gradient [K m-1]
T(z,t) = temperature at a certain depth, z, during a certain time, t [oC]
T(z) = average annual temperature at a certain depth, z [oC]
T0 = average surface temperature [
o C]
A0 = Amplitude of temperature at the surface [
o C]
Az = Amplitude of temperature at depth, z [
o C]
D = damping depth [m]
𝜔𝑦 = yearly angular frequency [rad s-1]
t = time [s]
Ф = phase angle [rad]
λ = thermal conductivity [W m-1 K-1]
C = volumetric heat capacity [J m-3 K-1]
ρ = density [kg m-3]
c = specific heat [J kg-1 K-1]
𝑧𝑎𝑠𝑙 = depth of the under limit of active surface layer [m]
𝑧𝑝𝑏 = depth of the permafrost base [m]
𝑇𝑚𝑎𝑥 = maximum monthly average temperature [oC]
𝑇𝑚𝑖𝑛 = minimum monthly average temperature [oC]
T1956 = average temperature at the year 1956 [
o C]
Tyr = average temperature at a certain year [
o C]
yr = year [-]
dTrad = temperature difference caused by radiation [
o C]
Toff = temperature offset [
o C]
Tfluc = temperature fluctuation [
o C] β = slope angle [o] φ = latitude [oN] ψ = aspect [o] elev. = elevation [m]
Page 7 of 57 2.3 Depth permafrost base
Figure 1: the permafrost layers and its temperature gradient in a homogenous soil.
In Figure 1 the temperature is schematic schematically illustrated, it shows that the annual mean surface temperature, T0, can be obtained by extrapolating the geothermal gradient (Osterkamp,
2003). The thermal gradient is shown as the seasonally invariant temperature gradient in Figure 1. Likewise, the base of the permafrost can be obtained by extrapolating the thermal gradient from the annual mean surface temperature, T0, until the depth where the temperature reaches 0.0 °C. In the
absence of a thermal gradient, the temperature at a certain depth and time can be obtained by the following function (van Boxel, 1986):
𝑇(𝑧, 𝑡) = 𝑇0+ 𝐴0 exp (− 𝑧
𝐷) sin (𝜔𝑦∗ 𝑡 + Ф − 𝑧
𝐷) [m] (1)
Different variables and their meaning are provided in the nomenclature. In this paper, the thermal gradient is introduced with its depth and integrated into this function. The thermal gradient comes forward out of the law of Fourier:
x = 𝛥𝑇 𝛥𝑧 = 𝐵 −𝜆 [K m -1 ] (2) Therefore the function becomes:
𝑇(𝑧, 𝑡) = 𝑇0+ 𝑥 ∗ 𝑧 + 𝐴0 exp (− 𝑧
𝐷) sin (𝜔𝑦∗ 𝑡 + Ф − 𝑧
𝐷) [m] (3)
For an average temperature in depth sin (𝜔𝑦𝑡 + Ф − 𝑧
𝐷) = 0, implying the following equation:
𝑇(z) = 𝑇0+ 𝑥𝑧. As the temperature at the base of the permafrost is equal to zero, 𝑇(z) = 0, the
function can be rewritten: 𝑧𝑝𝑏= −
𝑇0
𝑥 [m] (4)
2.4 Homogenous soil
Equation (4) is perfectly suitable for a homogenous soil. This means that the thermal gradient is the same for the whole soil. However, there is no such thing as a homogenous soil. The pedology –
Page 8 of 57
the pedogenesis, morphology and soil classification – is very complex and diverse. Although this complexity, the thermal gradient differs only from 0.025 to 0.030 K m-1 in the top of the earth crust (Osterkamp & Romanovsky, 1999). Therefore, in this paper the assumption was made that the thermal gradient, 𝑥, is equal to 0.0275 K m-1
throughout the soil.
Besides, another assumption was made, that the thermal conductivity, 𝜆, is set at 2.59 W m-1 K-1. As the thermal conductivity of clay is 2.92 W m-1 K-1, 𝜆 of ice is 2.10 W m-1
K-1 (De Vries & van Wijk, 1963). This paper made the assumption that every soil is a homogenous soil, that consists of 60% of clay minerals and 40% of ice. Whenever 𝜆 is divided between clay and ice, λ, is equal to 2.92*0.4 + 2.10*0.6 = 2.59 W m-1 K-1 (Monteith, 1973; Verkerk et al., 2008). The thermal conductivity of clay minerals lays in between the values of 𝜆 of other soil components, therefore it was chosen. The thermal conductivity of other soil components are shown in Table 1 (De Vries & Wijk, 1963).
Soil component Thermal conductivity [W m-1 K-1]
Quartz 8.80 Clay minerals 2.92 Organic matter 0.25 Water 0.57 Ice 2.10 Air (20 oC) 0.025
Table 2: the thermal conductivity of soil components (De Vries & van Wijk, 1963).
2.5 Depth active surface layer
Figure 1 displays that the temperature fluctuation becomes less in depth (van Boxel, 1986), at a certain depth the temperature amplitude approaches zero (van Boxel, 1986). In this paper, the depth where the amplitude is equal or less to 0.1 oC is seen as the upper limit of seasonally invariant temperature. Therefore, the depth of this point is seen as the depth of the active surface layer. The amplitude of temperature at a certain depth is equal to the 𝐴0 exp (−
𝑧
𝐷) (van Boxel, 1986), out of the
equation (1). Thus, the amplitude at a certain depth is approached with the follow equation: 𝐴z =
𝐴0 𝑒𝑥𝑝 (− 𝑧
𝐷). If Az ≤ 0.1
o
C then the temperature fluctuations are seen as not significant, making it that: 𝑧𝑎𝑠𝑙= −𝐷 ∗ ln ( 0.1 𝐴0 ) [m] (5) with 𝐷 = √( 2∗𝜆 𝐶∗𝜔𝑦) [m] (6) 𝐶 = 𝑓𝑟𝑎𝑐𝑡𝑖𝑜𝑛𝑐𝑙𝑎𝑦∗ 𝜌𝑐𝑙𝑎𝑦∗ 𝑐𝑐𝑙𝑎𝑦+ 𝑓𝑟𝑎𝑐𝑡𝑖𝑜𝑛𝑖𝑐𝑒∗ 𝜌𝑖𝑐𝑒∗ 𝑐𝑖𝑐𝑒 [J m-3 K-1] (7) 𝜔𝑦= 2𝜋 60∗60∗24∗365≈ 1.99E-7 [rad s -1 ] (8) Additionally, the amplitude of temperature at the surface, A0, is assumed to be equal to the amplitude
of monthly average temperature:
A0 =
𝑇𝑚𝑎𝑥−𝑇𝑚𝑖𝑛
2 [
o C] (9) As the paper assumed that the soil was a homogeneous one consisting of 60% clay and 40% ice. The volumetric heat capacity is equal to:
Page 9 of 57 𝐶 = 0.6 ∗ 2.7𝐸3 ∗ 0.9𝐸3 + 0.4 ∗ 0.9𝐸3 ∗ 2.2𝐸3 ≈ 2.2𝐸6 J kg-1
K-1 (Monteith, 1973; Verkerk et al., 2008) making it that the damping depth, 𝐷 ≈ √( 2∗2.59
3.41∗1.99𝐸−7) ≈ 3.4 m.
2.6 Depth permafrost table
The depth of the permafrost table represents the thickness of the active layer, this is also shown in Figure 1. Therefore, is the permafrost table the lower limit of the active layer.
3 Methods
3.1 Descriptive
The study area encompasses the Swiss Alps. This paper made use of the dataset of PermaNet (http://www.alpine-permafrostdata.eu/evidences-inventory). This dataset contains an inventory of permafrost evidence, that shows the locations/site points and conditions where the existence of the permafrost is examined. The ID, site names, coordinated, vegetation, elevation, slope angle and aspect of site points are withdrawn from the PermaNet dataset. What is more, only site points that have the occurrence of permafrost labeled as ‘yes’ and ‘unknown’ were included. If a site point had missing values, with the exception of vegetation, this site point was excluded.
The current report focusses on mudflow, landslides, and slumps hazards rather than on weathering. Therefore, the slope angle of ≤ 45o
was chosen because mudflows mainly consist of sand, clay and some debris (O’Brien & Julien, 1988). And saturated sand has the highest angle of repose that is equal to 45o (Rahardjo et al., 1995). At an angle greater than 45o it is likely that the whole soil just contains hard bedrock (WRB, 2006). Therefore, only site points were included with a slope angle ≤ 45o
.
3.2 Predicting air temperature
For 51 weather stations in the Swiss Alps elevation, latitude, longitude, and climatological values (1941-1970) were obtained for the mean annual temperature, temperature of the warmest month and temperature of the coldest month were obtained from Wernstedt (1972). All withdrawn data was converted to SI units: meter, decimal degrees, and oC.
A multiple regression model was built to predict the average air temperature. Prespecified variables were entered in the model using forward selection with entry criterium set at 0.05. The prespecified variables were elevation, latitude, and longitude. This analysis was performed in MATLAB R2017b software. The threshold for statistical significance was set at a nominal p = 0.05 level.
Since temperature data from Wernstedt (1972) is from 1941-1970, the outcome of regression predicts the average air temperature of 1956. Air temperature has risen by 0.1 oC per decade since the late 1950s (IPCC, 2013), making it that:
Tyr = T1956 + 0.010 * (yr - 1956) [
o
C] (10) For every site point the average air temperature of 2018 was calculated.
Page 10 of 57 3.3 Average surface temperature
To calculate the active surface layer and the active layer, the air temperature of 2018 needed to be converted into the average surface temperature:
T0 = Tyr + dTrad [
o
C] (11) with,
dTrad = Toff + Tfluc* (cos(β) cos(φ) + …
sin(β) sin(φ) cos(π - ψ)) [oC] (12)
The equation (10) and (11) were withdrawn from van Boxel (2002), the Toff is equal to the net
radiation at the ground converted into temperature. The net radiation at the ground is estimated to be 35 W*m-2 around in the Alps region, this estimation is based Raschke (2005) that used data of the International Satellite Cloud Climatology Project from 1991 to 1995. The net radiation at the ground was converted into temperature, to make the temperature offset out of the equation (12): Toff =
35*0.0021=0.0735 oC (Šafanda, 1999). The Tfluc = 1.5
o
C was chosen so that it is most in line with the radiation difference of a slope and aspect described in Šafanda (1999). Afterward, for every site point the average surface temperature of 2018 was calculated.
3.4 Temperature amplitude
Succeeding, temperature amplitude was calculated of the withdrawn Wernstedt (1972) data. The average monthly air temperature between 1941-1970 was given. The annual air temperature amplitude was calculated by equation (9). Then, again a multiple regression model was built to predict the air temperature amplitude for every site point. Hereby, two assumptions were made:
the air temperature amplitude does not change over time, and the amplitude of temperature at the surface is equal to the amplitude of air temperature. Prespecified variables were entered in the model using forward selection with entry criteria set at 0.05. The model with a statistical significance level of p< 0.05 was used.
3.5 Calculations
Thereafter, four main calculations were made in excel based on the described equations in the framework.
3.5.1 Depth of active layer
To calculate the depth of the permafrost table and therefore the depth of the active layer a derivative of equation (1) was used. The depth of the permafrost table can only be calculated during the summer when the average annual maximum temperature is reached. Therefore, sin (𝜔𝑡 + Ф −
𝑧
𝐷) = 1 (van Boxel, 1986), the whole active layer is thawed, and the temperature of the permafrost
table is equal to 0.0 oC. This makes the equation to: 𝑇0+ 𝑥𝑧 + 𝐴0 exp (−
𝑧
𝐷) = 0 [
o
C] (13) This can be rewritten as:
𝑧𝑝𝑡= − 𝐷 ∗ ln (−
𝑇0+𝑥∗𝑧𝑝𝑡
A0 ) [m] (14)
In the equation (14) 𝑧𝑝𝑏 is the estimate and the solution. Therefore, iterative calculation is needed. If
the solution is correct, the estimate of 𝑧𝑝𝑏 going into the calculation has to be equal to the answer 𝑧𝑝𝑡
Page 11 of 57
until the answer and the estimate are the same. Figure 2 illustrates how this is done in Excel. The depth of the permafrost table is found after three iterations. Figure 2 also shows whenever permafrost is absent, the depth of the permafrost is equal to zero after 3 iterations. 10 iterative calculations were done to be certain. Moreover, for the first estimation of 𝑧𝑝𝑡, Z0 in Figure 2, x is neglected as it is low
x=0.0275 K m-1. Therefore the equation (14) becomes solvable:
𝑧𝑝𝑡= −𝐷 ∗ ln (−
𝑇0
𝐴0 ) [m] (15)
This equation (15) was used to have a reasonable first estimate.
Figure 2, the iteration number on the x-axis and the estimated permafrost table depth [m] on the y-axis, given the permafrost table depth by iterative calculation.
3.5.2 Depth active surface layer
Equation (5) was used to calculate the depth of the active surface layer.
3.5.3 Depth permafrost base
Equation (4) was used to calculate the depth of the permafrost table.
3.5.4 Alteration of temperature
In order to discuss how the distribution of permafrost changes in Switzerland, as the air temperature will increase with 2.0 oC according to the climate medium emission scenario by 2100 (Gobiet et al., 2014). Therefore the average annual surface temperature, T0, was increased by 2.0
o C for all calculations concerning the depth of the permafrost base and the permafrost table. To provide a result what will happen the thickness of the permafrost and its active layer. So it is assumed that the amplitude of temperature fluctuations will not change by 2100, but T0 will increase.
Page 12 of 57
4 Results
4.1 Descriptive
62 site points are encountered of the PermaNet evidence table. The permafrost occurrence is labelled as ‘unknown’ for 38 site points and as ‘yes’ for 24 site points. The coordinates of all site points are ranging from 45.9 oN - 46.8oN in the latitude, and from 10.0 oE – 7.1oE in the longitude. 2% of the
vegetation has complete coverage, 16% partly coverage, 10% sparse, 60% none, and 12% has a missing value. The elevation is raging from 1580 – 4450 m (mean µ 2647 m, standard deviation (SD) 388.8m). The median slope is 27.5o (mean µ 27.2o, SD 11.4o). The aspect has a mean of 183.5o and a median of 180o.
4.2 Predicting average surface temperature
The temperature data was available for 51 Swiss weather stations. Table 3 provides a summary of the observed Wernstedt (1972) data.
Variables Range Mean Median Standard deviation
Elevation [m] 237.00 – 2499.00 953.70 685.00 612.80 Latitude [oN] 45.83 – 47.70 46.83 46.82 0.46 Longitude [oE] 6.15 – 10.30 8.35 8.50 0.99 Average air temperature [oC] -2.28 – 11.78 6.36 7.22 3.31 Temperature amplitude [oC] 7.00 – 10.31 8.75 8.89 0.72 Temperature amplitude without outliers [oC] 7.00 – 9.72 8.65 8.89 0.66
Table 3: a summary of the observed Wernstedt (1972) data in suitable units.
Estimated coefficients Estimate SE p-val
β0 [°C] 66.8 9.1 2.0E-9
β1 [°C/m] -5.1E-3 1.4E-4 1.1E-37
β2 [°C/°N] -1.2 0.2 1.7E-7
DF: 51 Root mean squared error: 0.59 Adjusted R-Squared: 0.97 F-statistic: 752 p-value : 5.93E-37
Table 4: multiple regression output of used model, temperature_1956 = β0 + β1elevation + β2latitude.
The average air temperature is best predicted by the elevation and latitude:
T1956 = 66.8 – 5.1E-3*elev. – 1.2*φ [
o
C] (16) Table 4 provides the results of the predictors used in the model. In addition, Figure 3 shows that the multiple regression with elevation and latitude predicts the average air temperature better than a linear regression model with only elevation as predictor. In the linear regression, 10 values are outside the +/- 1 oC confidence interval, this decrease to only 5 values with the multiple regression. Making it that the adjusted R-squared is higher with multiple regression than the R-squared with linear regression. The air temperature is modelled for every site point, converted into air temperature of 2018 (10) and surface temperature of 2018 by the formula (11). This modelled surface temperature, T0, is
Page 13 of 57
Figure 3, modelled temperature versus observed temperature. The drawn black line is the 1:1 line, the dashed lines shown a difference of +1. 0 °C (red) or -1.0 °C (blue) between modelled temperature and observed temperature.
4.3 Temperature amplitude
The temperature amplitude, A0, at every weather station was calculated and is shown in Figure 4 (blue and red dots) and the values are given in Appendix 2.
The air temperature amplitude can be predicted by:
A0 = 9.4 – 7.2E-4*elev. [
o
C] (17) The predictor variables latitude and longitude are non-significant to the model.
As Figure 4 shows, some calculated A0 are outliers (red dots), confirmed by the boxplot.
Outliers are defined as values that are lower than the 25th percentile minus 1.5 multiplied by interquartile range (Q1 – 1.5*IQ), and values that are higher than the 75th percentile minus 1.5
multiplied by the interquartile range (Q3 + 1.5*IQ) (Ghasemi & Zahediasl, 2012). IQ is defined as the 75th percentile minus the 25th percentile. By the removal of the outliers the Root Mean Squared Error decreases from 0.576 to 0.339 and the R-squared increases from 0.375 to 0.740. Therefore the model without outliers is used to calculate A0, the output of this model is shown in Table 5.
Page 14 of 57
Figure 4, temperature amplitude, A0, versus elevation, with the boxplot of the residuals for the yellow line. Yellow line itself
shows the regression lines including outliers. While the Green line shows the regression lines excluding outliers. The boxplot illustrates that there are four outliers for the residuals. The red line is the median (-0.108 oC), blue boxplot shows the interquartile range (Q3= -0.319 oC,
Q3 = 0.168 oC), within the whiskers there are no outliers (upper limit = 0.412 oC and lower limit= -0.562 oC).
Estimated coefficients Estimate SE p-value
β0 9.5 8.9E-2 8.8E-56
β1 -9.2E-4 1.4E-4 9.4E-15
DF: 47 Root mean squared error: 0.34 R-squared: 0. 74 F-statistic: 128 p-value : 9.4E-15
Table 5, output regression temperature_amplitude_without_outliers = 9.5 -9.2E4*elevation
4.4 Calculations
4.4.1 Depth active layer
The depth of the permafrost table, and therefore the depth of the active layer, was calculated for 2018 and 2100 for every site point (Figure 5). In 2018 41 out of the 62 site points had permafrost. Between 2330 m and 2660 m elevation, the occurrence of permafrost depended on other variables, like slope angle, aspect or A0. Above the elevation of 4450 m the site point had permafrost, but without an active layer, because the permafrost table was calculated at 0 m. In 2100 only 18 out of the 62 site points had permafrost. The active layer increased in thickness for every site point that still contained permafrost, except for the site point at 4450 m elevation. The occurrence of permafrost became dependent from 2536 m till 3100 m elevation.
Page 15 of 57
Figure 5, shows the depth of the permafrost table in 2018 (left side) and in 2100 (right side). On the x-axis is the elevation in meters and on the y-axis is the depth of permafrost table in meters.
4.4.2 Depth active surface layer
Figure 6 shows the depth of the active surface layer. The zasl was dependent of the A0. As A0 is
dependent on the elevation, zasl was also dependent on elev., which made a correlation of 1.000. with
the formula:
zasl= -3.2E-8*elev.
2
-2.7E-4*elev. +1.547E1 [m] (18)
As assumed that the temperature amplitude does not change over time, the active surface layer did not change over time either. R-squared is 0.954 for a linear regression between elevation elev. and zasl.
Figure 6, the depth of the active surface layer, with a 1.000 correlation between elevation (x-axis) in meters and active surface layer (y-axis) in meters. This is the active surface layer in 2018 and 2100.
4.4.3 Depth permafrost base
The depth of the permafrost base was calculated for the year 2018 and 2100. As given before, in 2018 41 site points had permafrost and in 2100 only 18 site points had permafrost left. There seemed to be a correlation between the depth of permafrost base and elevation. The depth of the permafrost will decrease of time making it that, in 2100 the depth of the permafrost will be less than in 2018.
y = -3,2E-08x2 - 2,7E-04x + 1,5E+01
R² = 1,000 13,4 13,6 13,8 14 14,2 14,4 14,6 14,8 15 15,2 1500 2000 2500 3000 3500 4000 4500 z, act iv e su rf ace laye r [m ] Elevation [m]
Page 16 of 57
Figure 7, the depth of permafrost of 41 site points in 2018 and of 18 site point in 2100, with a trend line and the correlation formula.
5 Discussion
Climate change is bringing along all kind of hazards, like mass movement originating from permafrost in the susceptible dense Alps. This paper provides an insight on how climate change affects the distribution of permafrost, the thickness of the active surface layer and the active layer in
Switzerland. It contributes to the lack of knowledge the current probability models have about the permafrost and its active layer, and how it is projected and evolved in line with current climate change model (Boeckli et al., 2002). Therefore it could be a base for the hazard policymakers.
This study used data from PermaNet to collect the characteristics of 62 different site points throughout the Swiss Alps and collected the air temperature from 52 weather stations. Based on the data of the weather stations, mathematical models were developed to predict the air temperature and temperature amplitude for the 62 site points. Thereafter the depth of permafrost, active surface layer and the active layer of the 62 site points was calculated for the current situation (2018) and the situation in 2100.
According to these calculations, permafrost will vanish for 56% of all site points in 2100. No permafrost was found below 2330 m in 2018, and below 2536 m in 2100. Between 2330-2600 m the occurrence of permafrost depended on other variables, in 2100 this ranged will shift to 2536 to 2660 m according to the calculations. At the site points where permafrost still occurs in 2100, the depth of the active layer will increase compared to 2018. For site points that will still have permafrost in 2100, the permafrost layer will become thinner compared to 2018, as the permafrost base decreases in its depth and the active layer increases in depth. Given that the surface temperature amplitude was assumed not to change over time, the influence of climate change on the active surface layer could not be well predicted.
Figure 5 shows that in 2018, permafrost is found from the elevation of 2330 m onward. Between the range of 2330 to 2600 m the occurrence of permafrost is depended on other variables. This range will shift to a range of 2536 to 2660 m in 2100 according to the calculations, meaning that no permafrost is found beneath 2536 m elevation. Therefore permafrost will vanish in 56% of all site points that had permafrost. This is quite in line with Guo & Wang (2017) who predict 47% worldwide
Page 17 of 57
decrease of permafrost within 100 years with the same climate scenario. Lunardini (1996), Gavrilova (1993), and Nelson & Anisimov (1993) also predict that permafrost will vanish in the next 100 years at many places.
Figure 5 also illustrates that for every site that still contains permafrost and has an active layer, the active layer will grow with time. Likewise, Haeberli & Gruber (2009) stated that the active layer will increase. The active surface layer, however, stays the same over time illustrated in Figure 6. This is explained by the fact that A0 is the only variable that influences the depth of the active surface
layer. Given that A0 was assumed constant over time, the depth of the active surface layer did not
change either. Meanwhile, the depth of the permafrost decreased as Figure 7 shows. The latter was also seen on the Qinghai-Tibet plateau (Wu et al., 2010).
This whole paper is based on Figure 1 that shows the thermal gradient through a homogenous soil, shown as the red lines. Based on the results this red line will move to the right because of climate change. This is supported by Changwei et al. (2015), who monitored the ground temperature to 30 m depth. When the average surface temperature increased, the temperature at a certain depth increased as well, making it that temperature line (red line in Figure 1) shifts to the right.
Although, the results of the current paper are in line with the results of others some limitation should be mentioned. First of all, the reader should bear in mind that the study is based on a model and like every model it is a simplification of reality (Brodsky & Watson 2000). In this model, the assumption was made that the soil is homogeneous until the base of the permafrost, causing the thermal gradient to be constant with depth. In addition, this paper assumed that the soil contained of 60% clay particles and 40% ice for every site point, while the pedology and geology of the Alps are very complex (Meyre et al., 1999). Therefore, it is unlikely that the thermal gradient is the same for every site point (Boeckli et al., 2012).
Secondly, the thermal conductivity is also based on 60% clay particles and 40% ice. However, this chosen thermal conductivity is low compared to the conductivity of hard bedrock (Khan, 2002), as shown in Table 6. Hard bedrock is often found within 25 cm of the surface in mountain areas (WRB, 2006). Therefore, hard bedrock and its higher conductivity could not be neglected. With a higher conductivity (Khan, 2002), the damping depth (6) will increase, making the active surface layer (5) and the active layer (14) to become thicker, and the thermal gradient to become less (2) resulting in a lower depth of the permafrost base (4). Nevertheless, in this current paper the thermal gradient is not only underestimated by the latter, it is also overestimated. As a change in thermal gradient is seen just above the permafrost base (Osterkamp & Burn, 2003), this is caused by the difference in the thermal conductivities between the frozen and unfrozen ground (Osterkamp & Burn, 2003). An explanation can be found in the fact that clay has a 39% higher thermal conductivity compared to ice (Table 2), therefore the thermal gradient just above the permafrost base is 39% less (2). The damping depth decreases and the depth of the active surface layer, the depth of the active layer, and the depth of permafrost base with it. However, the effect of this change in thermal gradient is not seen in the measurements of Osterkamp & Romanovsky (1999) or in the monitoring of Changwei (2015).
Page 18 of 57 Type of concrete Thermal conductivity (W m-1 K-1)
Basalt 4.03 – 4.30 Limestone 3.15 – 3.49 Siltstone 3.52 – 5.22 Quartzite 8.58 – 8.63 Dolerite 2.43 Haemite 5.44
Table 6: the thermal conductivity of concrete types (Khan, 2002).
To calculate the temperature of the site points, regression models for the air temperature and temperature amplitude were built. In these models elevation, latitude, and longitude were the
prespecified variables. The multiple regression had an adjusted R-squared of 0.97 for air temperature and an R-squared of 0.74 of temperature amplitude, meaning that there is still some unexplained variance for both regression models, and thus other variables besides the prespecified variables might explain the variance better.
To convert the obtained air temperature to surface temperature, some assumptions were made. As the Wernstedt (1972) data gives the average monthly air temperature between 1941-1970, it is assumed that the average of air temperature presents a value that is representative for 1956. This air temperature is converted to 2018 (10) and then converted into surface temperature (11).
Toff is assumed to be a constant variable that is used to calculate the dTrad (12). The net
radiation at the ground was estimated to be 35 W m-2 around in the Alps region. However, it is within a range of 24 to 72 W m-2 (Raschke et al., 2005). The offset of 35 W m-2 was chosen, because it is more in line with colder regions with a lower net radiation at the ground than warmer regions with a higher net radiation at the ground. Otherwise, 48 W m-2 would have been a more obvious choice because it is the middle of the range 24 to 72 W m-2. What is more, the cloud effect of -12 to -24 W m -2
described by Raschke (2005) is not taken in consideration. The effects of cloud, however, do influence the net radiation at the ground (Raschke et al., 2005). The effects of clouds might change in time as well under the influence of climate change. In addition concerning the Toff, the net radiation at
the ground was converted into temperature by multiplying it with 0.0021.
The Tfluc is chosen to be most in line with Šafanda (1999). Šafanda (1999) provides a more
precise estimation of chosen Tfluc of 1.5
o
C, as the temperature difference between 0o to 50o in slope is 1.37 oC, and the temperature difference between 00 to180o in aspect is 1.40 oC (Šafanda, 1999). Furthermore, the thesis does not engage factors such as soil moisture content and biomass of vegetation cover (Du et al., 2016), that also influence the surface temperature. In addition, surface temperature strongly depends on snow conditions rather than the atmospheric temperatures alone, but snow is a hardly predictable interface (Haeberli & Gruber, 2009).
What is more, this paper assumed that the temperature amplitude at the surface is equal to the air temperature amplitude. While, surface temperature and it amplitude are also dependent on
variables such as radiation and precipitation (Boeckli, 2012).
This all causes that for 10 out of the 24 there is no permafrost found while the evidence table of PermaNet said otherwise.
In addition, it was assumed that the amplitude of temperature will be the same in 2100. However, Qian et al., (2011) found that in China the amplitude of the surface temperature
significantly (P < 0.01) decreased by -4.6% or -0.63 oC between 1961-2007. A decrease is also found for the Alps region by Gobiet et al. (2014), who predict the spatially averaged change is seasonally varying between -20.4% in summer and +10.4% in winter until the end of the 21st century, resulting in
Page 19 of 57
a decrease of the annual temperature amplitude. This would mean that the active surface layer would be different in 2100.
A side note has to be made to Figure 6, the R-squared=1.000 is provided for the equation (18), to show that zasl can be calculated by elev.. However, Spiess and Neumeyer (2010) claim that the
evaluation of R-squared is an inadequate measure for nonlinear models. Therefore the R-squared is calculated for a linear regression between elev. and zasl (trendline shown in Appendix 1). The
R-squared is calculated for a linear regression between elev. and zasl is equal to 0.954, making it that
there is a very high correlation between elev. and zasl.
Despite that the paper contains several assumptions, the trends found are in line with the existing literature and therefore it is a good base for further research. Further research should focus on determining the pedology and geology of every site point or area to make a more precise estimation of the active surface layer, active layer and the depth of the permafrost. Moreover, further research needs to examine more closely the links between the variable radiation, precipitation, and snow cover to predict A0 and T0 and how those latter two change over time. This all might be useful for a permafrost
probability map in the Alps, or might be encountered into hazard maps for mudflows, landslides and rockfalls.
6 Conclusion
The current research suggests that the average surface temperature is best predicted by elevation and latitude and that amplitude of temperature at the surface is best predicted by elevation only. According to the performed calculations, permafrost will vanish in a homogeneous soil in many places. In 2100, when the air temperature will be increased by 2 °C, there will not be any permafrost below 2536 m elevation. Nowadays, the occurrence of permafrost between 2330 m and 2660 m elevation depends on other variables, like slope angle or aspect, according to the performed
calculations. In 2100, this range will shift to between 2536 m and 3100 m elevation. At the site points where permafrost still occurs in 2100, the depth of the active layer will increase compared to 2018. For site points that will still have permafrost in 2100, the permafrost layer will become thinner compared to 2018, as the permafrost base rises to lower depths and the active layer increases in depth. Given that the surface temperature amplitude was assumed not to change over time, the influence of climate change on the active surface layer could not be well predicted.
In short, the red line, displayed in Figure 1, that represents the temperature at a certain depth and time, will shift to the right over time. When the red line shifts to the right, permafrost will vanish at some places and the depth of the active layer will increase at points where permafrost still occurs.
7 References
Anisimov, O. A., Shiklomanov, N. I., & Nelson, F. E. (1997). Global warming and active-layer thickness: results from transient general circulation models. Global and Planetary Change, 15(3-4), 61-77.
Boeckli, L., Brenning, A., Gruber, S., & Noetzli, J. (2012). Permafrost distribution in the European Alps: calculation and evaluation of an index map and summary statistics. The Cryosphere, 6(4), 807
Brodsky, D., & Watson, B. (2000, May). Model simplification through refinement. In Graphics Interface (Vol. 2000, pp. 221-228).
Page 20 of 57
Changwei, X., Gough, W. A., Lin, Z., Tonghua, W., & Wenhui, L. (2015). Temperature-dependent adjustments of the permafrost thermal profiles on the Qinghai-Tibet Plateau, China. Arctic, antarctic, and alpine research, 47(4), 719-728.
Crozier, M. J. (1999). Prediction of rainfall‐triggered landslides: A test of the antecedent water status model. Earth Surface Processes and Landforms: The Journal of the British Geomorphological Research Group, 24(9), 825-833.
De Vries, D. A. (1963). Thermal properties of soils. pp 210–235. In WR van Wijk (ed.) Physics of plant environment. North-Holland Publishing, Amsterdam.
Du, E., Zhao, L., Wu, T., Li, R., Yue, G., Wu, X., ... & Wang, Z. (2016). The relationship between the ground surface layer permittivity and active-layer thawing depth in a Qinghai–Tibetan Plateau permafrost area. Cold Regions Science and Technology, 126, 55-60.
Gavrilova, M. K. (1993). Climate and permafrost. Permafrost and Periglacial Processes, 4(2), 99-111. Ghasemi, A., & Zahediasl, S. (2012). Normality tests for statistical analysis: a guide for
non-statisticians. International journal of endocrinology and metabolism, 10(2), 486.
Gobiet, A., Kotlarski, S., Beniston, M., Heinrich, G., Rajczak, J., & Stoffel, M. (2014). 21st century climate change in the European Alps—a review. Science of the Total Environment, 493, 1138-1151.
Guo, D., & Wang, H. (2017). Permafrost degradation and associated ground settlement estimation under 2 C global warming. Climate Dynamics, 49(7-8), 2569-2583.
Gruber, S., & Haeberli, W. (2007). Permafrost in steep bedrock slopes and its temperature‐related destabilization following climate change. Journal of Geophysical Research: Earth Surface, 112(F2).
Gruber, S., Hoelzle, M., & Haeberli, W. (2004). Rock‐wall temperatures in the Alps: modelling their topographic distribution and regional differences. Permafrost and Periglacial Processes, 15(3), 299-307.
Haeberli, W. (2005). Climate change and glacial/periglacial geomorphodynamics in the Alps: a challenge of historical dimensions. Geografia Fisica e Dinamica Quaternaria, Suppl, 2, 9-14. Haeberli, W., & Gruber, S. (2009). Global warming and mountain permafrost. In Permafrost soils (pp.
205-218). Springer, Berlin, Heidelberg.
IPCC (2013) Summary for Policymakers. In: Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change [Stocker, T.F., D. Qin, G.-K. Plattner, M. Tignor, S.K. Allen, J. Boschung, A. Nauels, Y. Xia, V. Bex and P.M. Midgley (eds.)]. Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA.
Kääb, A. (2008). Remote sensing of permafrost‐related problems and hazards. Permafrost and Periglacial Processes, 19(2), 107-136.
Khan, M. I. (2002). Factors affecting the thermal properties of concrete and applicability of its prediction models. Building and Environment, 37(6), 607-614.
Lunardini, V. J. (1996). Climatic warming and the degradation of warm permafrost. Permafrost and Periglacial Processes, 7(4), 311-320.
Meyre, C., De Capitani, C., Zack, T., & Frey, M. (1999). Petrology of high-pressure metapelites from the Adula nappe (Central Alps, Switzerland). Journal of Petrology, 40(1), 199-213.
Monteith, J. L. (1973). Principles of Environmental Physics. London: Edward Arnold (publishers) Limited.
Nelson, F. E., & Anisimov, O. A. (1993). Permafrost zonation in Russia under anthropogenic climatic change. Permafrost and Periglacial Processes, 4(2), 137-148.
Page 21 of 57
O'Brien, J. S., & Julien, P. Y. (1988). Laboratory analysis of mudflow properties. Journal of hydraulic engineering, 114(8), 877-887.
Osterkamp, T. E., & Romanovsky, V. E. (1999). Evidence for warming and thawing of discontinuous permafrost in Alaska. Permafrost and periglacial Processes, 10(1), 17-37.
Osterkamp, T.E.; Burn, C.R. (2003) "Permafrost", in North, Gerald R.; Pyle, John A.; Zhang, Fuqing, Encyclopedia of Atmospheric Sciences, 4, Elsevier, pp. 1717–1729, ISBN 0123822262 Qian , C., Fu, C., & Wu, Z. (2011). Changes in the amplitude of the temperature annual cycle in
China and their implication for climate change research. Journal of Climate, 24(20), 5292-5302.
Rahardjo, H., Lim, T. T., Chang, M. F., & Fredlund, D. G. (1995). Shear-strength characteristics of a residual soil. Canadian Geotechnical Journal, 32(1), 60-77.
Raschke, E., Ohmura, A., Rossow, W. B., Carlson, B. E., Zhang, Y. C., Stubenrauch, C., ... & Wild, M. (2005). Cloud effects on the radiation budget based on ISCCP data (1991 to 1995).
International Journal of Climatology: A Journal of the Royal Meteorological Society, 25(8),
1103-1125.
Šafanda, J. (1999). Ground surface temperature as a function of slope angle and slope orientation and its effect on the subsurface temperature field. Tectonophysics, 306(3-4), 367-375.
Scherler, M., Hauck, C., Hoelzle, M., Stähli, M., & Völksch, I. (2010). Meltwater infiltration into the frozen active layer at an alpine permafrost site. Permafrost and Periglacial Processes, 21(4), 325-334.
Spiess, A. N., & Neumeyer, N. (2010). An evaluation of R 2 as an inadequate measure for nonlinear models in pharmacological and biochemical research: a Monte Carlo approach. BMC
pharmacology, 10(1), 6.
Su, Z. (2002). The Surface Energy Balance System (SEBS) for estimation of turbulent heat fluxes. Hydrology and earth system sciences, 6(1), 85-100.
van Boxel, J. H. (2002). Modelling Gloabal Radiation for the Portofino area in Italy. Amsterdam: Institute for Biodiversty and Ecosystem Dynamics (IBED) University of Amsterdam.
van Boxel, J. H. (1986). Heat Balance Investigations in Tidal Areas (Doctoral dissertation, Centre for Mathematics and Computer Science).
Verkerk, G., Broens, J. B., de Groot, P. A. M., Kranendonk, W., Vogelezang, M. J., Westra, J. J., Wevers-Prijs, I. M., (2008). BINAS havo/vwo. Groningen/Houten: Noordhoff Uitgevers bv. Wernstedt, F. L. (1972). World Climatic Data. Climatic Data Press, Lemont, Pennsylvania, United
States, 523 pp.
WRB 2006. World reference base for soil resources – a framework for international classification, correlation and communication. World Soil Resources Reports 103, International Union of Soil Sciences, ISRIC – World Soil Information and Food and Agriculture Organization of the United Nations, Rome.
Wu, Q., Zhang, T., & Liu, Y. (2010). Permafrost temperatures and thickness on the Qinghai-Tibet Plateau. Global and Planetary Change, 72(1-2), 32-38.
Page 22 of 57
Appendix
1 extra graphs
Figure 1: active surface layer versus elevation with a linear trendline.
Figure 2: depth permafrost table versus elevation in 2100 seems to have a rational trendline: zpt = 1800 /(elev.-2400).
2 spreadsheet of data and calculations the spreadsheet can be found at:
https://docs.google.com/spreadsheets/d/1rQiPnFyBjYGT3XUuz3Nlb-5uIcLDJ6CnwcFdvhfRm2o/edit?usp=sharing
3 matlab script
%% The influence of climate change on permafrost and its
active layer in the Swiss Alps.
%% J.E.D. Kerkmeijer
%% thesis MATLAB script %%
clear
all
Page 23 of 57
%% 1. regression test mean air temperature as predictor
and elevation,
% longitude, and latitude as prespecified variables
(Wernstedt, 1972)
Elevation= [406
456
463
1441
1864
275
237
572
2072
1711
1331
1141
633
1561
959
1018
433
670
405
480
2475
450
568
446
1077
990
705
685
553
1380
239
276
498
605
412
1090
1377
1331
1786
1985
2499
508
474
Page 24 of 57
1252
549
470
702
2095
912
1609
515];
Elevation_range=max(Elevation)-min(Elevation);
Elevation_mean=mean(Elevation);
Elevation_median=median(Elevation);
Elevation_SD=std(Elevation);
Latitude=[47.38
46.88
47.38
46.63
46.78
47.63
46.18
46.95
46.50
46.55
46.30
47.03
46.85
46.80
46.92
46.82
47.57
46.80
46.20
47.03
45.83
47.70
46.68
47.65
46.98
47.10
47.35
46.78
46.52
46.38
46.17
46.02
47.05
Page 25 of 57
46.72
46.43
47.03
46.67
46.47
47.05
46.37
47.25
47.05
46.90
46.80
46.23
47.22
47.42
46.55
47.28
46.03
47.38
];
Latitude_range=max(Latitude)-min(Latitude);
Latitude_mean=mean(Latitude);
Latitude_median=median(Latitude);
Latitude_SD=std(Latitude);
Longitude=[8.05
8.65
9.55
8.60
9.68
7.67
9.03
7.43
9.17
9.88
9.13
6.88
9.53
9.83
9.18
8.42
8.90
7.15
6.15
9.07
7.17
8.45
Page 26 of 57
7.85
9.18
6.60
6.83
7.77
7.78
6.63
7.63
8.80
8.97
8.30
8.20
6.92
8.78
8.85
8.25
8.50
6.98
9.33
9.43
8.25
10.30
7.35
7.53
9.38
8.57
8.92
7.75
8.53];
Longitude_range=max(Longitude)-min(Longitude);
Longitude_mean=mean(Longitude);
Longitude_median=median(Longitude);
Longitude_SD=std(Longitude);
Mean_air_temperature_1956=[8.28
9.11
8.39
3.00
2.72
9.11
11.50
8.11
0.50
1.39
6.39
5.61
Page 27 of 57
8.22
2.78
5.89
5.39
8.28
7.72
9.89
7.61
-1.72
8.39
8.11
8.39
4.72
6.00
6.72
7.22
9.22
5.22
11.78
11.28
8.78
8.00
10.00
5.22
4.78
4.00
2.28
1.78
-2.28
8.78
8.39
4.89
9.78
8.50
6.89
-0.22
7.11
3.61
8.78];
Mean_air_temperature_1956_range=max(Mean_air_temperature_
1956)-min(Mean_air_temperature_1956);
Mean_air_temperature_1956_mean=mean(Mean_air_temperature_
1956);
Mean_air_temperature_1956_median=median(Mean_air_temperat
ure_1956);
Page 28 of 57
Mean_air_temperature_1956_SD=std(Mean_air_temperature_195
6);
%linear model between Elevation and air temperature:
m1=fitlm(Elevation,Mean_air_temperature_1956)
% Linear regression model:
% y ~ 1 + x1
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% __________ __________ _______
__________
%
% (Intercept) 11.359 0.20439 55.576
6.3214e-46
% x1 -0.0052431 0.00018081 -28.998
1.6512e-32
%
%
% Number of observations: 51, Error degrees of freedom:
49
% Root Mean Squared Error: 0.784
% R-squared: 0.945, Adjusted R-Squared 0.944
% F-statistic vs. constant model: 841, p-value = 1.65e-32
%multiple regression model Elevation and Latitude
predictors
m2=fitlm([Elevation,Latitude],Mean_air_temperature_1956)
% Linear regression model:
% y ~ 1 + x1 + x2
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% __________ _________ _______
__________
%
% (Intercept) 66.806 9.0654 7.3694
2.0009e-09
% x1 -0.0055074 0.0001436 -38.353
1.1444e-37
% x2 -1.1787 0.19268 -6.1172
1.6598e-07
%
%
Page 29 of 57
% Number of observations: 51, Error degrees of freedom:
48
% Root Mean Squared Error: 0.593
% R-squared: 0.969, Adjusted R-Squared 0.968
% F-statistic vs. constant model: 752, p-value = 5.93e-37
%multiple regression model Elevation, Latitude and
Longitude as predictors
m3=fitlm([Elevation,Latitude,Longitude],Mean_air_temperat
ure_1956)
% Linear regression model:
% y ~ 1 + x1 + x2 + x3 + x4
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% __________ __________ _______
__________
%
% (Intercept) 67.932 9.2944 7.3089
3.1419e-09
% x1 -0.0056082 0.00017583 -31.896
5.0957e-33
% x2 -1.2181 0.20087 -6.0642
2.3204e-07
% x3 0.071277 0.089454 0.7968
0.42966
% x4 0.00016336 0.00020431 0.7996
0.42805
%
%
% Number of observations: 51, Error degrees of freedom:
46
% Root Mean Squared Error: 0.598
% R-squared: 0.97, Adjusted R-Squared 0.967
% F-statistic vs. constant model: 371, p-value = 2.36e-34
%this model has two non-significant estimated
coefficients (predictors).
%% 2. regression test: temperature amplitude as predictor
and elevation, longitude, latitude and precipitation as
prespecified variables
Temperature_amplitude=[8.89
8.50
9.31
8.64
7.67
Page 30 of 57
8.92
9.72
9.19
7.75
10.31
7.89
8.00
8.94
9.31
8.86
8.39
9.25
9.00
9.00
9.25
7.69
9.28
8.50
9.11
8.36
8.44
8.75
9.00
8.86
8.33
8.94
9.56
9.06
9.11
8.64
8.50
8.11
9.78
7.19
7.19
7.00
8.89
9.11
10.11
9.67
9.33
8.89
7.58
8.25
9.06
8.92];
Page 31 of 57
Temperature_amplitude_range=max(Temperature_amplitude)-min(Temperature_amplitude);
Temperature_amplitude_mean=mean(Temperature_amplitude);
Temperature_amplitude_median=median(Temperature_amplitude
);
Temperature_amplitude_SD=std(Temperature_amplitude);
%linear regression between Elevation and Temperature
amplitude
m4=fitlm(Elevation,Temperature_amplitude)
% Linear regression model:
% y ~ 1 + x1
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% ___________ __________ _______
__________
%
% (Intercept) 9.4324 0.15032 62.748
1.7907e-48
% x1 -0.00072059 0.00013298 -5.4188
1.8137e-06
%
%
% Number of observations: 51, Error degrees of freedom:
49
% Root Mean Squared Error: 0.576
% R-squared: 0.375, Adjusted R-Squared 0.362
% F-statistic vs. constant model: 29.4, p-value =
1.81e-06
%so to see how it looks in an scatterplot:
%equation:
%y=9.4324-0.00072059*Elevation;
x=[0:1:2499];
y=9.4324-0.00072059*x;
figure(2)
scatter(Elevation,Temperature_amplitude,
'b'
,
'filled'
);
xlabel(
'Elevation [m]'
);ylabel(
'Temperature amplitude
[oC]'
); title(
'Elevation vs. Temperature amplitude'
)
hold
on
plot(x,y,
'g'
)
Page 32 of 57
Residuals=Temperature_amplitude-(9.4324-0.00072059*Elevation);
figure(3)
boxplot(Residuals);ylabel(
'difference modelled and
observed A0 [oC]'
);xlabel(
'amount of boxplots'
)
;title(
'boxplot residuals: temperature amplitude = 9.432
– 7.206E – 4*elevation'
)
%median
median(Residuals)
%interquartile range
Residuals_IQR=iqr(Residuals)
%upperlimit boxplot
Residuals_Q1=quantile(Residuals,0.75)
%lowerlimit boxplot
Residuals_Q3=quantile(Residuals,0.25)
%whiskers
whiskerdown=Residuals_Q1-1.5*Residuals_IQR
whiskerup=Residuals_Q3+1.5*Residuals_IQR
%provides 4 outliers:
Temperature_amplitude_outliers=[ 10.31
10.11
9.78
9.31];
Elevation_outliers=[ 1711
1252
1331
1561];
%outliers excluded:
Elevation_without_outliers= [406
456
463
1441
1864
275
237
572
2072
1331
1141
633
959
1018
433
670
405
Page 33 of 57
480
2475
450
568
446
1077
990
705
685
553
1380
239
276
498
605
412
1090
1377
1786
1985
2499
508
474
549
470
702
2095
912
1609
515];
Temperature_amplitude_without_outliers=[8.89
8.50
9.31
8.64
7.67
8.92
9.72
9.19
7.75
7.89
8.00
8.94
8.86
8.39
9.25
9.00
Page 34 of 57
9.00
9.25
7.69
9.28
8.50
9.11
8.36
8.44
8.75
9.00
8.86
8.33
8.94
9.56
9.06
9.11
8.64
8.50
8.11
7.19
7.19
7.00
8.89
9.11
9.67
9.33
8.89
7.58
8.25
9.06
8.92];
Temperature_amplitude_without_outliers_range=max(Temperat
ure_amplitude_without_outliers)-min(Temperature_amplitude_without_outliers);
Temperature_amplitude_without_outliers_mean=mean(Temperat
ure_amplitude_without_outliers);
Temperature_amplitude_without_outliers_median=median(Temp
erature_amplitude_without_outliers);
Temperature_amplitude_without_outliers_SD=std(Temperature
_amplitude_without_outliers);
%linear regression Elevation, Air temperature without
outliers
m5=fitlm(Elevation_without_outliers,Temperature_amplitude
_without_outliers)
Page 35 of 57
% y ~ 1 + x1
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% ___________ _________ _______
__________
%
% (Intercept) 9.483 0.088772 106.82
8.7571e-56
% x1 -0.00091649 8.099e-05 -11.316
9.4034e-15
%
%
% Number of observations: 47, Error degrees of freedom:
45
% Root Mean Squared Error: 0.339
% R-squared: 0.74, Adjusted R-Squared 0.734
% F-statistic vs. constant model: 128, p-value = 9.4e-15
%equation without outliers:
temperature_amplitude=9.483-0.00091649*Elevation;
y2=9.483-0.00091649*x;
%Figure that provides scatter plot indicating outliers
and regression lines
%including and excluding the outliers:
figure(4);
scatter(Elevation_without_outliers,Temperature_amplitude_
without_outliers,
'b'
,
'filled'
)
hold
on
scatter(Elevation_outliers,Temperature_amplitude_outliers
,
'r'
,
'filled'
)
plot(x,y,
'k'
)
plot(x,y2,
'g'
)
xlabel(
'Elevation [m]'
);ylabel(
'Temperature amplitude
[oC]'
);title(
'scatterplot A0 versus
elevation'
);legend(
'Temperature amplitude'
,
'temperature
amplitude outliers'
,
'regression with
outliers'
,
'regression without outliers'
)
hold
off
Longitude_without_outliers=[8.05
8.65
9.55
8.60
9.68
Page 36 of 57
7.67
9.03
7.43
9.17
9.13
6.88
9.53
9.18
8.42
8.90
7.15
6.15
9.07
7.17
8.45
7.85
9.18
6.60
6.83
7.77
7.78
6.63
7.63
8.80
8.97
8.30
8.20
6.92
8.78
8.85
8.50
6.98
9.33
9.43
8.25
7.35
7.53
9.38
8.57
8.92
7.75
8.53];
Latitude_without_outliers=[47.38
46.88
47.38
Page 37 of 57
46.63
46.78
47.63
46.18
46.95
46.50
46.30
47.03
46.85
46.92
46.82
47.57
46.80
46.20
47.03
45.83
47.70
46.68
47.65
46.98
47.10
47.35
46.78
46.52
46.38
46.17
46.02
47.05
46.72
46.43
47.03
46.67
47.05
46.37
47.25
47.05
46.90
46.23
47.22
47.42
46.55
47.28
46.03
47.38];
Page 38 of 57
%multiple regression without outliers with Elevation and
Longtitude as
%preditors:
m6=fitlm([Elevation_without_outliers,Longitude_without_ou
tliers],Temperature_amplitude_without_outliers)
% Linear regression model:
% y ~ 1 + x1 + x2
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% ___________ __________ _______
__________
%
% (Intercept) 9.2918 0.45311 20.507
3.8015e-24
% x1 -0.00091803 8.1812e-05 -11.221
1.6993e-14
% x2 0.023365 0.054273 0.43051
0.66893
%
%
% Number of observations: 47, Error degrees of freedom:
44
% Root Mean Squared Error: 0.342
% R-squared: 0.741, Adjusted R-Squared 0.729
% F-statistic vs. constant model: 63, p-value = 1.23e-13
%multiple regression without outliers with Elevation and
Latitude as
%predictors
m7=fitlm([Elevation_without_outliers,Latitude_without_out
liers],Temperature_amplitude_without_outliers)
% Linear regression model:
% y ~ 1 + x1 + x2
%
% Estimated Coefficients:
% Estimate SE tStat
pValue
% ___________ __________ _______
__________
%
% (Intercept) 18.843 5.0721 3.715
0.00057005
% x1 -0.00095976 8.2318e-05 -11.659
4.7736e-15
Page 39 of 57
% x2 -0.19897 0.10781 -1.8456
0.071688
%
%
% Number of observations: 47, Error degrees of freedom:
44
% Root Mean Squared Error: 0.33
% R-squared: 0.759, Adjusted R-Squared 0.748
% F-statistic vs. constant model: 69.2, p-value =
2.62e-14
%% 3. Data-analysis, withdrawn data from PermaNet
evidence table
elevation=[2489
2410
2430
2519
2135
2860
3000
2700
2394
2501
2805
2370
2760
2820
2780
2420
2580
2575
2900
2570
2690
3100
2360
3295
2477
2474
2516
2455
2330
4450
2600
2100
2250
Page 40 of 57
2400
2500
2500
2020
2941
2410
2946
2680
2410
2765
2690
2740
2630
2500
2549
2558
2538
2536
1580
2820
2840
2394
2660
3077
3410
3410
2910
2910
2670];
elevation_min=min(elevation);
elevation_max=max(elevation);
elevation_mean=mean(elevation);
elevation_SD=std(elevation);
latitude=[46.13675018
46.15222292
46.15222292
46.13659868
46.19472189
46.49583737
46.41000114
46.42086605
46.7479251
46.74687972
46.49593327
46.77716126
Page 41 of 57
