Spatial modelling of climatic water balance using GIS methods

Agnieszka Wypych

Spatial modelling of climatic water balance using GIS methods

MSc thesis under the supervision of dr. Katarzyna Ostapowicz

MSc thesis submitted in the framework of, and according to the requirements of the UNIGIS Master of Science programme (Geographical Information Science & Systems)

Jagiellonian University, Kraków, Paris Lodron University of Salzburg

Krakow 2012

2

I declare that all sources used in the thesis were properly acknowledged.

The thesis is fully my work and it was not and will not be submitted as a thesis elsewhere.

……… ………

Date Signature

3

ABSTRACT

Recent developments in GIS techniques have produced a wide range of powerful methods for capturing, modelling, and displaying of climate data. The main aim of this study is to find the best spatialisation method to describe spatial differentiation of climatic water balance (CWB) in Poland. Monthly mean values of air temperature and precipitation totals for 60 synoptic stations as well as monthly totals of solar radiation (21 measuring points) were taken into consideration. Source material covered period 1985-2006.

Regarding the prior research as well as data availability the potential evapotranspiration data was calculated by Turc formula.

CWB was conducted with two methods simultaneously: simple and multiply linear regression (with latitude, altitude and distance from the coast line as variables) and map algebra. Validation showed map algebra as the best spatialisation method. Nevertheless the obtained results proved also that except for the method local factors are of the great importance in CWB modelling especially in the mountains and at the coast. To optimize the method it is necessary to reduce the research scale using more in-situ data what would enable to include more local variables as land form and land cover into the analyses.

Key words:

spatial analysis, regression model, map algebra, climatic water balance, Poland

4

T

1. Introduction ……… 5

1.1. GIS in climatology and meteorology……… 5

1.2. Climatic water balance estimation – methods overview .… 8 2. Aims ……… 11

3. Study area ……… 12

4. Data ……… 14

5. Methods ……… 18

5.1. Spatialisation methods ……… 20

5.1.1. Regression models: simple linear regression, multiple linear regression ……… 22 5.1.2. Map algebra ……… 26

5.2. Validation ……… 33

6. Spatial differentiation of the climatic water balance in Poland …… 35

6.1. Assessment of climatic water balance spatial differentiation in Poland using regression models ……… 37

6.2. Assessment of climatic water balance spatial differentiation in Poland using map algebra ……… 42

6.3. Validation of CWB models for Poland ……… 47

7. Discussion ……… 56

8. Conclusions ……… 61

9. References ……… 63

10. List of figures ……… 64

11. List of tables ……… 66

5

1. I

Recent development of GIS techniques has produced a wide range of powerful methods for modelling and displaying of climatological and meteorological data.

Advanced processing methods allow for detailed analysis of climate elements at different temporal and spatial scales. GIS techniques designed to temperature and atmospheric precipitation fields have received the most attention so far. However, researchers are often interested not in meteorological elements themselves but in the information that can be extracted from them in a form of various indices, which are useful in environmental and human sciences (Tveito et al. 2008).

1.1. GIS in climatology and meteorology

In the year 2001 within the frame of the European Union and European Cooperation in the Field of Scientific and Technical Research (COST) Action “The use of Geographic Information Systems in Climatology and Meteorology” was established (COST 719). The main objective of the Action was to broaden and enhance the potential of GIS in the fields of climatology and meteorology by developing applications in those fields, with emphasis on the procedures and capabilities for integrating and adding value to data from various sources and on quality control and presentation of climate and other related data (Dobesch et al. 2001). Within 5 years after the end of the project many methods which were worked out are still widely used in many detailed analysis. Most of the present studies concerning using GIS in climatology and meteorology are results of the project either directly or experience-based at least (Chapman, Thornes 2003, Dyras et al.

2005).

All the studies concerning GIS in climatology and meteorology can be divided into three main groups of research:

1) investigation of data standards, interoperability, and accessibility,

2) implementation of interpolation methods of different meteorological elements, 3) GIS applications (Chapman, Thornes 2003, Dyras et al. 2005).

6 Within the first group also some studies concerning an estimation of contemporary use and farther possibilities of GIS in climatology and meteorology can be included (Dobesch et al. 2001, Christakos et al. 2002, Tveito, Schöner 2002, Chapman, Thornes 2003, Wel et al. 2004, Almarza, Gutierrez 2004, Madej 2004, Hoffman 2004, Dyras et al.

2004a; Wilhelmi 2004, Woolf et al. 2005, Ustrnul 2006, Dobesch et al. 2007). All the authors confirm that although GIS are widely accepted within the atmospheric science but currently their role is still limited mainly in meteorology.

Within the second group conducted studies give the complex information about spatialisation techniques used so far in climatology mapping (Dobesch et al. 2001, Tveito, Schöner 2002, Stach, Tamulewicz 2003, Dobesch et al. 2007). The most valuable publications are by H. Dobesch, O. E. Tveito and P. Bessemoulin (2001) as well as O. E. Tveito and W. Schöner (2002). They have presented and described all (or almost all) interpolation techniques used in climate and weather research. According to the authors the role of GIS in meteorology and climatology is incontestable however it is necessary to assimilate particular applications for the atmospheric science research (Dahlström 2001, Dobesch et al. 2001, Tveito, Ustrnul 2003, Dobesch et al. 2007).

The last group consists of particular applications where GIS were used to solve detailed research problems. In most publications spatial differentiation of thermal conditions (Álvaro Pimpăo 2003, Ustrnul, Czekierda 2003, 2005, Ustrnul 2006) or precipitation was presented (Sobik et al. 2001, Sobolewski 2001, Bac-Bronowicz 2002, 2003, Loukas et al. 2004). Also snow cover (Almarza, Gutierrez 2004), drought extension (Szalai, Bihari 2004, Dolinar, Sušnik 2004) as well as floods possibility (Todini et al.

2004) were the topics. GIS tools appeared to be useful also in forest climatology determine macroclimate conditions and describe spatial variability and distribution of microclimates (Durło 2004) and agroclimatology where the areas of potential frosts around the wine- growings were detected (Madelin, Beltrando 2005). GIS applications were also widely used to gain and process satellite and radar data (Dyras et al. 2004b, Cros et al. 2004, Dyras, Serafin-Rek 2005).

Those application-specific results have shown that – what was previously mentioned – researchers are often interested not in climatological or meteorological elements themselves but in the information that can be extracted from them (Tveito et al.

2008).

7 One of such indices, presenting a basis of a climatic assessment of water resources of a given area, is a climatic water balance (CWB). It focuses on the difference between the precipitation (P) and evapotranspiration (E). The knowledge of spatial distribution of the climatic water balance seems to be very important according to its comprehensive application in spatial management, agriculture and hydroclimatological modelling, nevertheless it has not been a subject of detailed analyses so far, mainly because the climatic water balance and its spatialisation is very complex. It is, the first and foremost, a subject associated with a problem of spatial interpolation of evapotranspiration, which varies considerably with changes in the natural environment, secondly, it relates to a problem of data availability. Given the complicated nature of the subject, it is no wonder that there exist many GIS methods that attempt to model the spatial differentiation of evapotranspiration (e.g. Nováky 2002, Xinfa et al. 2002, Kar, Verma 2005, Loheide, Gorelick 2005, Fernandes et al. 2007, Vicente-Serrano et al. 2007). Also remote sensing techniques are becoming more commonly used to address this research issue and often supplements ground-based observations (Woolhizer, Wallace 1984, Rosema 1990, Kalma et al. 2008).

As abovementioned, the climatic water balance is the complex index dependant on many different variables such as solar radiation, relief, land use, urban development and others. As the data availability is poor – especially regarding evapotranspiration – the problem of appropriate index interpretation mainly due to its spatial differentiation arises.

GIS techniques and methods give possibilities to enclose different data processing and integration methods, with complex analyses and modelling methods.

8 1.2. Climatic water balance estimation – methods overview

The climatic water budget was introduced in the middle of 20^th century by Thornthwaite (1948). He described the budget as the balance of precipitation, potential evapotranspiration, and actual evapotranspiration, taking into account both soil moisture utilization and soil moisture recharge (Oliver, Fairbridge 1987). According to Thornthwaite and his colleagues (Thornthwaite, Mather 1957), an average climatic water budget model can be expressed using two interrelated equations:

P = ET + S, (1)

and

PE = ET + D, (2)

where P – precipitation, ET – evapotranspiration, PE – potential evapotranspiration, S – moisture surplus, D – moisture deficit. The first equation defines water inflow, outflow, and storage; the second equation describes energy demands. The procedure designed by Thornthwaite and Mather (1957) to calculate climatic water balance parameters is still widely used in CWB research (e.g. Tateishi, Ahn 1996, Kar, Verma 2005, Vicente-Serrano et al. 2010).

An evapotranspiration process is the principal component of the climatic water balance, as it returns 60% to 80% of precipitation back into the atmosphere. In order to determine the value of the CWB index, the magnitude of evapotranspiration must be properly estimated. Owing to the difficulty of obtaining accurate field measurements, evapotranspiration is commonly computed from weather data using empirically derived formulas (Jaworski 2004). A large number of more or less empirical methods and procedures have been developed over the last 50 years and are designed to estimate evaporation, potential evaporation and evapotranspiration from different climatic variables (Jaworski 2004). One of the most cited is the Thornthweite’s and Holtzmann’s formula defining evaporation from lakes using air density, specific humidity and wind speed. To describe potential evapotranspiration as the function of sunshine duration and temperature another Thornthweite’s formula is widely used mainly in the United States whereas Turc’s procedures, with EP calculated using solar radiation and temperature, are often used and

9 mostly well known in Europe (Jaworski 2004). Some of the methods are only valid under specific climatic and agronomic conditions and cannot be applied under conditions different from those under which they were originally developed (Allan et al. 2004).

Currently the FAO Penman-Monteith method (3) is recommended as the standard method for the definition and computation of the reference evapotranspiration (ETo).The reference evapotranspiration provides a standard to which evapotranspiration at different periods of the year or in different regions can be compared (Allan et al. 2004).

^{( )}

( )

( ) (3)

where:

ET_o – reference evapotranspiration [mm day^-1], R_n – net radiation at the crop surface [MJ m^-2 day^-1], G – soil heat flux density [MJ m^-2 day^-1],

T – mean daily air temperature at 2 m height [°C], u₂ – wind speed at 2 m height [m s^-1],

e_s – saturation vapour pressure [kPa], ea – actual vapour pressure [kPa],

es - ea – saturation vapour pressure deficit [kPa],

 – slope vapour pressure curve [kPa °C^-1],

 – psychrometric constant [kPa °C^-1].

Contemporary climate change has brought detailed analysis of climate elements mainly the ones responsible for droughts as major causes of agricultural, economic and environmental damage. Thereby the detailed studies of precipitation and evapotranspiration process become essential. To detect, monitor and analyse droughts a new drought index: the Standardised Precipitation-Evapotranspiration Index (SPEI) which is nothing more that the standardized difference between the elements has been introduced (Vicente-Serrano et al. 2010). The SPEI allows to compare drought severity through time and space, since it can be calculated over a wide range of climates (The Standardised Precipitation...).

In Poland, most evapotranspiration and climatic water balance research is focused on the identification of a model that would best suit weather conditions in Poland. All the

10 investigation have shown that evapotranspiration process play a decisively role in the water cycle as well as energy exchange at the active surface (Jaworski 2004). At the same time they emphasize the lack of homogenous data, rare distributed measuring points as well as great complexity of evapotranspiration process what make it more difficult to create any model. However, some physically based mathematical models, with less of more complicated formulas, have been developed. They can be implemented under different physiographic, vegetation and climate conditions (Jaworski 2004). The following formulas were used in existing research: Turc (1964) – potential evapotranspiration, Bac (1970) reference evaporation – local index, and Penman modified to Polish conditions (Sarnacka et al. 1983). All three methods were applied to the analysis of the measurements data (Wild scale and GGI-3000 pan evaporimeter). Although, the Turc's method produced the largest differences between evapotranspiration totals measured in situ and values derived empirically (also shown by Nováky 2002 for Hungary), the method proved to be useful because of data availability issues. It was also selected because of other research that has shown that it is best at determining relationships with elevation (Kowanetz 1998, Kowanetz 2000).

Since 2007 the Drought Monitoring System for Poland has been provided by the Institute of Soil Science and Plant Cultivation – State Research Institute in Puławy (Institute of Soil…). In the system, meteorological conditions that are causing drought are evaluated by the climatic water balance expressed by the difference between the precipitation and potential evapotranspiration (by Penman formulae). Nevertheless it has not been the subject of detailed analysis so far.

11

2. A

The aim of the study is to develop new methodology for climatic water balance index implementation with use of geographical information systems (GIS) in cases when there is no appropriate spatial information given from in-situ observations.

Using GIS enables to involve different data processing, its integration from different sources and complex synthesis. The application to analyse such a complex index as climatic water balance (CWB), seems to be well justified and leads to acquire the best, i.e. regionally and locally specified spatial information.

My specific objectives are: (1) to find out which spatialisation method is the most appropriate to be used regarding both complex index of climatic water balance and also the failure of the data, (2) to assess an error magnitude as well as advantages and defects of all the methods implemented, (3) to state which explanatory variables (predictors) could be used according to the research area and time (with monthly and annual resolution) scale.

All seem to be really important when long-term analyses are taken into account as well as climate modelling.

12

3. S

The study area covers whole territories of Poland which is located in the Central Europe (Fig. 3.1). The country has an area of 313 000 km2 and occupies the eastern part of the North European Plain between the Baltic Sea in the north and the Carpathian Mountains and the Sudeten Mountains in the south. Poland’s relief has a parallel latitudinal scheme of the principal types of landscape (Fig. 3.1). The average elevation in Poland is 173 m a.s.l. with the lowest point situated in the Vistula River Delta (-1.8 m a.s.l.) and the highest in the Tatra Mountains (Rysy; 2499 m a.s.l.). About 91% of the country’s territory lies below 300 m a.s.l., which indicates a typical lowland-type landscape. Poland’s location and its relief features are important in the formation of its weather and climate patterns.

The beltlike layout of the main landscape regions is favourable to free zone circulation. For this reason, one can observe the collision of oceanic and continental air masses in this area (Ustrnul, Czekierda 2009). At the same time, the southern part of the country possesses relatively diverse relief, which impacts local weather and climate conditions. This also appears on the Polish coast where the Baltic Sea activity is apparent in many weather events. The factors mentioned above as well as the lie of the land conditions make it possible to observe significant weather differences across Polish territory. Given the diverse weather conditions present in Poland, one can state that Poland lies in a moderately warm transitional climate zone that is the product of colliding air masses over its territory.

Poland’s climate is characterized by considerable variability in weather conditions and significant seasonal differences from year to year. This can be easily inferred from the character of winter seasons, which can be of an oceanic or continental type. The influence of local conditions, especially relief and atmospheric circulation, is expressed in the variability of climate elements, mainly precipitation (Woś 2011).

13 Fig. 3.1. Location of Poland in Europe

The lie of the land as well as country location induce the opinion that the analysis based on the testing area (Poland) seems to be entirely justified and quite representative for the broader region e.g. the central-eastern Europe.

14

4. D

Although climatic water balance is quite a simple index to compute, complexity of evapotranspiration process due to responsible factors makes it really difficult to obtain exact values of CWB for meteorological (current) analyses. However for climatological purposes (long term or general analyses) the value of ETP with a high precision can also be evaluated by the use of simplified models, which include meteorological elements that are normally measured by the meteorological stations.

Therefore analyses of the climatic water balance (CWB) are usually developed for regions where input data, mainly air temperature and precipitation, can be readily obtained from meteorological stations. The research presented herein is based on mean monthly values of air temperature and precipitation totals obtained from 60 meteorological stations (Fig. 4.1A, 4.1B) as well as monthly totals of solar radiation coming from 21 solar stations (Fig. 4.1C). The data have been covered respectively for periods 1951-2006 and 1986- 2006.

Fig. 4.1. Location of the source data stations: A) temperature, B) precipitation, C) solar radiation

Not all solar stations collect necessary meteorological data therefore for the detailed analyses of climatic water balance only the data from 16 stations and covering the period from 1986 to 2005 were used (Fig. 4.2).

15 Fig. 4.2. Location of the climatic water balance data sources stations

Meteorological data were completed with information about the elevation and land cover.

To generate relief information the Shuttle Radar Topography Mission Digital Elevation Model (SRTM DEM) – the result of the international project spearheaded by the National Geospatial-Intelligence Agency (NGA) and the National Aeronautics and Space Administration (NASA) – was used (Earth Resources Observation…). The data is available as 3 arc second (approximately 90 meters resolution) DEMs with the vertical error of, as being reported, less than 16 meters (Earth Resources Observation…). The model was resampled into 90 meters spatial resolution and transformed to PUWG-92 coordinate system and finally stored in *.img data format (Fig. 4.3).

CORINE Land Cover (CLC) database was used to obtain the land use information (European Environment Agency…). The database contains both raster and vector maps (with the minimum mapping unit of 25 ha) of land cover and land use made as the visual interpretation of satellite data (Landsat and SPOT). The distinguished polygons are grouped into 44 land use classes, 31 of which can be found in Poland. In the current study

16 the data generated for the year 2000 (CLC2000) were used. Four land use classes were distinguished: artificial surfaces (aggregated CLC level 1), agricultural areas (aggregated CLC level 2), forest and semi natural areas (aggregated CLC level 3), water bodies and wetlands (aggregated CLC level 4 and 5) (Fig. 4.4). The CLC map was transformed from original coordinate system the ETRS 1989 LAEA L52 M10 to PUWG-92 coordinate system and stored in *.shp data format (Fig. 4.4).

Fig. 4.3. Shuttle Radar Topography Mission Digital Elevation Model for Poland (SRTM DEM)

Fig. 4.4. Land use in Poland based on CLC 2000

17 Due to effective data processing for the detailed analysis both the SRTM DEM and the CLC 2000 were resampled into 1 km spatial resolution. The resampling methods used were those widely recommended for these kind of dataset. Nearest neighbour method was implemented for the CLC 2000 as it is used primarily for discrete data, such as a land use classification. The SRTM DEM was resampled with cubic convolution (bicubic) method, which is appropriate for continuous data. Interpolation tends to smooth out peaks and valleys, however the bicubic method creates quite an accurate reproduction of the statistical behaviour of the local terrain variability retaining local extremes to a decent degree (ArcGIS Desktop Help…).

18

5.

Although the climatic water balance (CWB) is a complex index dependant on many variables it was created based on a simplified definition where it is interpreted as the difference between the precipitation total (P) and potential evapotranspiration (PE). Given the limited nature of the source data Turc’s formula (1964) was used to obtain potential evapotranspiration values. This method is confirmed (Kowanetz 1998, 2000) to be suitable for describing the relationships between evapotranspiration and relief. The final formula was the following:

(4)

where:

P – monthly precipitation totals [mm]

t – monthly average air temperature [°C]

I – monthly sum of total solar radiation [cal cm^-2 day^-1]

Climatic water balance modelling was carried out using two approaches simultaneously: linear regression method and map algebra. Results of both were validated with common error estimators. The CWB values calculated for the 16 stations were used as the reference data (Fig. 4.2, Tab. 5.1).

The analysis was conducted for growing season treated as the period from May until October which is one of the most frequent used in agrometeorology.

19 Table 5.1. Geographical parameters of selected reference stations ( – longitude;  – latitude; H – altitude; d – distance from the coast, CLC – land use)

reference stations

  H D CLC

º

º

Kołobrzeg 15.58 54.18 3 0 112 / artificial surfaces Łeba 17.53 54.75 2 1 231 / agricultural areas Piła 16.75 53.13 72 133 121 / artificial surfaces Toruń 18.60 53.04 69 147 141 / artificial surfaces Mikołajki 21.58 53.78 127 140 313 / forest and semi

natural areas Koło 18.66 52.20 115 239 211 / agricultural areas Warszawa 20.96 52.16 106 266 124 / artificial surfaces Legnica 16.21 51.19 122 327 211 / agricultural areas Sulejów 19.87 51.35 188 336 211 / agricultural areas Jelenia Góra 15.79 50.90 342 349 112 / artificial surfaces Śnieżka 15.74 50.74 1603 366 324 / forest and semi

natural areas Kłodzko 16.61 50.44 356 416 211 / agricultural areas Bielsko-Biała 19.00 49.81 398 504 112 / artificial surfaces Zakopane 19.96 49.29 855 565 112 / artificial surfaces Kasprowy

Wierch 19.98 49.23 1991 572 332 / forest and semi natural areas Lesko 22.34 49.47 420 581 211 / agricultural areas

20 5.1. Spatialisation methods

The term interpolation is usually used alternatively to not strictly defined idea of spatialisation. It is, according to Dobesch et al (2007), a set of methods describing the dependency on neighbouring data of a dataset in a typically the Cartesian coordinates system. The main aim of the methods is to condense and visualize data however when used also to transform and derive the new data or extended information spatialisation becomes a part of spatial analysis.

All of the large variety of spatialisation methods are based on theoretical consideration, assumptions and constraints which must be fulfilled in order to use the method properly. According to Ustrnul and Czekierda (2003) the most adequate classification of spatial analysis methods should regards the one based on both theoretical and empirical fundamentals. Following the suggestion, methods of spatial analysis can be distinguished as deterministic, stochastic and combined ones that consist of both deterministic and stochastic procedures. Regardless the interpolation method, an assumption that a measured value gives the information about the place and, with less accuracy, its surroundings can be found at base (Magnuszewski 1999).

Deterministic methods are characterized by wide range of mathematical methods using both the easiest algorithms and also complex calculation procedures (Tveito, Schöner 2002). They use empirical models to explain spatial phenomena, so in fact are conceptual and less abstract than stochastic methods (Ustrnul, Czekierda 2005). Deterministic methods can be represented by the use of data from the nearest station (e.g. Thiessen polygons), the weighted linear combination of data from neighbouring stations (e.g.

Inverse Distance Weighting; IDW), approximate polynomial functions (trend surface analysis), exact polynomial functions, and radial basis functions (Dobesch et al 2007).

Stochastic methods usually apply probability theory and the concept of randomness in spatial processing. These methods allow the spatial co-variance to be included in the interpolation process. It means that interpolated surface is only the one of many that might have been observed, all of which could have produced the known data points (Tveito et al.

2008). Stochastic methods for spatial problems are often called as geostatistical methods.

Geostatistics theory is based on the idea that values measured at locations close together tend to be more alike than values measured at locations further apart. To quantify the spatial correlation which is then incorporated in the spatial model variogram or covariance

21 function are used. Geostatistical techniques have the capability of producing a prediction surface and also provide some measure of the certainty or accuracy of the predictions (Tveito, Schöner 2002).

The most known statistical method is kriging. It is based on both mathematical and statistical models and relies on the notion of autocorrelation. The rate at which the correlation decays can be expressed as a function of distance (ESRI 2003). It is assumed that the distance between points is related to their similarity what is presented on semivariogram (ESRI 2003).

Combined methods, what was previously mentioned, integrate the basis of stochastic and deterministic methods and are more and more often used in climatology and meteorology.

Besides abovementioned interpolation methods there are also many different where separate, new algorithms are being created or which are based on complex and very complicated mathematical and statistical interpolation assumptions.

In this study climatic water balance modelling was conducted with two methods simultaneously. Regarding correlations between environment elements linear regression was used in the first approach. Statistical relations between CWB and geographical factors:

latitude, altitude and distance from the coast line were taken into account (Wypych et al.

2010). The second approach was based on data modelling with map algebra procedure implementation.

22 5.1.1. Regression models: simple linear regression, multiple linear regression

As abovementioned, first approach postulated using regression models: simple linear regression (SLR) and multiple linear regression (MLR). Close relationships between climatic water balance and geographical factors became the basis (Wypych et al. 2010) with longitude, latitude, altitude as well as distance from the coast and land cover as explanatory variables.

The global regression model is expressed by:

∑ (5)

where:

yi –predicted values

xik – explanatory variables (predictors)

0 – model intercept

k – linear regression coefficients

i – regression residuals

Because of the limited number of samples and also the least (from all the analysed predictors) correlation coefficient between CWB values and longitude (Table 5.2) as well as small differentiation of land cover at measuring points (Table 5.1) it was finally decided to exclude those variables from regression models (6-8) and not use them in further analyses:

Z(s) = β0 + β2 H(s) + ε(s) (6)

Z(s) = β₀ + β₁ (s) + β2 H(s) + ε(s) (7)

Z(s) = β₀ + β₁ (s) + β2 H(s) + β₃ d(s) + ε(s) (8) where:

Z(s) – dependant variable φ(s) – latitude,

H(s) – altitude (m a.s.l.)

d(s) – distance from the coast (m) ε(s) – regression residuals

23 Table 5.2. Pearson’s correlation coefficient values (R) between CWB index and selected predictors; values statistically significant at =0.05 level bolded

predictor* May June July Aug. Sept. Oct. May-

Oct.

 -0.131 -0.216 -0.126 -0.234 -0.193 -0.370 -0.212

 0.283 0.318 0.113 0.385 0.355 0.599 0.345

H 0.864 0.853 0.841 0.864 0.862 0.911 0.878

d 0.314 0.348 0.166 0.382 0.340 0.552 0.357

*  – longitude;  – latitude; H – altitude; d – distance from the coast

On the basis of described linear dependences and using regression method CWB values were calculated for grid points delimited with the spatial resolution of 1 km.

Fig. 5.1. Workflow for first approach with the use of linear regression

 – latitude; H – altitude; d – distance from the coast

24 In the last step of creation of the climatic water balance spatial differentiation map (Fig. 5.1) data interpolation method with Radial Basis Functions (RBF) was used RBFs is the exact interpolation technique, it takes into account the general tendencies and also local variability. In RBFs the surface must go through each measured sample point but the technique minimizes its total curvature. The delimited values can either exceed or be lower than the real ones (Fig. 5.2). Although RBFs is the most useful interpolator in the case of many data points, it gives worse results and should not be implemented when the parameter’s distribution varies within short distances or when the data quality is poor and doubtful (ESRI 2001).

Fig. 5.2. Radial Basis Functions (RBF) – interpolation scheme (source: ESRI 2001).

ArcGIS Spatial Analyst creates radial basis function for every data point within delimited neighbourhood. Shape of the obtained surface varies together with the distance from the points. Parameter value for the point (Fig. 5.3) is a weighted average of the surfaces ø1, ø2 and ø₃ (Fig. 5.3).

modeled surface

measured data

25 Fig.5.3. Radial Basis Functions (RBF) for selected points – model scheme (source ESRI 2001)

Research conducted hitherto (Wypych, Ustrnul 2010) has confirmed the applicability of RBF method also for CWB index spatialisation.

Basis function is calculated for each

point separately

vertical section for y=5

To get value at point y=5 and x=7 mean value of Ø1, Ø2, Ø3is

calculated

26 5.1.2. Map algebra

The second approach was based on map algebra application. Model here is treated as sequence of different sub-functions leading to the simulation of final spatial differentiation image (Urbański 2008). The term “algebra” describes the way of data manipulating with defined mathematical formulas operating on raster matrix (Fig. 5.4).

Fig. 5.4. Simplified map algebra scheme (source: http://www.rockware.com)

Modelling here means a process of raster data transformation using GIS tools. Map algebra was allowed to build complex expressions and process them as a single command.

For the present study map algebra was used to get the final CWB spatial differentiation map. In the first step a set of maps presenting spatial distribution of climatic water balance index components: air temperature, precipitation totals and solar radiation was created using in-situ data (Fig. 5.5)

27 Fig. 5.5. Workflow for second approach with the use of map algebra

Component maps of climatic water balance index were constructed according to methodology worked out so far by international research teams dealing with GIS implementation in meteorology and climatology (Chapman, Thornes 2003, Dobesch et al.

2007). The method commonly used and suggested as probably the best is kriging (Chapman, Thornes 2003, Dobesch et al. 2007). Required value estimation at interpolation points consists in calculation of weighted moving average of points surrounding (Magnuszewski 1999). Such a procedure ensures to use different explanatory variables (predictors) to make the interpolated predicted estimation more reliable.

As first CWB component (Fig. 5.5), temperature spatial differentiation maps were created using residual kriging method (Ustrnul, Czekierda 2005). It is one of the combined

28 methods, including both deterministic approach (spatial differentiation of a variable is explained with empirical models) and the stochastic one (concept of randomness in spatial processing). First, a multivariate linear regression was applied, then the residuals from the regression model were calculated for each station, and finally the residuals were spatially interpolated by the ordinary kriging method. The final map is a result of both sub-layers summation (compare 9-11).

Z(s) = β₀ + β₁ x_i(s) + β₂ x₂(s) + …+ βp x_p + ε(s) (9) Z*(s) = β₀ + β₁ x_i(s) + β₂ x₂(s) + …+ βp x_p (10)

*(s) = Z(s) – Z*(s) (11)

where:

Z(s) – dependant variable Z*(s) – modelled value xi(s) – predictors

ε(s) – regression residuals

ε*(s) – estimated random model residuals

Several geographic parameters, including elevation, latitude, longitude, and distance to the Baltic coast (for stations located within 100 km), were used as predictor variables for air temperature. Multivariate linear regression was used to develop the initial statistical relationship between air temperature and terrain. Elevation played the most important role.

Precipitation totals were interpolated over the area of Poland with ordinary kriging method (Łupikasza et al. 2005, Łupikasza 2007). The prediction obtained is a linear combination of measured values with weights depending on the spatial correlation between the data. It assumes intrinsic stationarity and data isotropy (spatial correlation function should be homogenous in all directions) as well as the lack of any trend. The predictions are weighted linear combinations of the available data (Tveito et al. 2008). Linear coefficient are calculated under the condition of a uniformly unbiased predictor and under the constraint of minimal prediction error variance (kriging variance).Ordinary kriging is an exact interpolation method used when the mean value for the whole series stays unknown (compare equation 12) – the algorithm uses local means instead (Namysłowska- Wilczyńska 2006, ESRI 2003). The disadvantage of the method is that meteorological variables can rarely be considered as an intrinsic stationary random process. Using the data

29 alone it is no way to decide if the observed pattern is the result of spatial autocorrelation (among the errors ε(s) with the constant µ) or also the trend.

Z(s) = µ(s) + (s) (12)

where:

µ(s) – unknown constant ε(s) – random residual

The solar radiation surface has been obtained by double-acting. Firstly solar radiation field was generated with ordinary kriging method for the point data from 21 measuring stations (compare to chapter 4). Simultaneously but independently the application solar analyst ArcGIS was used. All the necessary information was implemented, i.e. sun hours, attitude for the site area, radiation parameters (diffuse factor and transmittivity) and topographical factors: slope, aspect and shaded relief counted from SRTM. Because of element sensitivity to local conditions (astronomical, geographical and meteorological parameters) numbers of different settings in solar analyst application were tested to achieve satisfactory final results. Mainly: diffuse factor and transmittivity were adjusted. Mean diffuse factor, i.e. fraction of global normal radiation flux that is diffused, ranges from 0 to 1 and should be set according to atmospheric conditions. Typical values are 0.2 for very clear sky conditions and 0.3 for generally clear sky conditions (ESRI 2003). For the area of Poland reaches around 0.5 was received. Transmittivity is treated as a property of the atmosphere and is the ratio of the energy received at the upper edge of the atmosphere to that reaching the Earths's surface by the shortest path (in the direction of the zenith), averaged over all wavelengths. Values range from 0 (no transmission) to 1 (complete transmission). Typically observed values are 0.6 or 0.7 for very clear sky conditions and 0.5 for only a generally clear sky (ESRI 2003). For Poland value of 0.4 could be assumed. Different combinations testing showed that the results closest (with the least error) to the real data are those gained with the default settings (respectively 0.3 and 0.5).

Created images were compared to each other (Fig. 5.6a, Fig. 5.6b). Real data map, because of limited number of measuring points (only 21 measuring stations across the country irregular distributed), is characterized by significant errors when compared to reference data. Furthermore because of the lack of reference points in particular regions it is difficult to estimate the model quality herein. Whereas solar analyst model, besides the

30 significant difference in the mountain areas (Fig. 5.6a, Fig. 5.6b) which are the results of using terrain (STRM) as a variable, gives less errors over the most area of Poland.

Therefore in final analyses in-situ data was used only as the reference data to select proper parameters and model estimation.

Fig. 5.6a. Spatial differentiation of solar radiation totals [cal·cm^-2] in Poland within vegetation period (May /V/–July /VII/): SA – solar analyst, OK – ordinary kriging

31 Fig. 5.6b. Spatial differentiation of solar radiation totals [cal·cm^-2] in Poland within vegetation period (August /VIII/–October /X/): SA – solar analyst, OK – ordinary kriging

32 All the layers created were used as an input data to potential evapotranspiration model by Turc's formula (air temperature together with solar radiation map) and further on with precipitation map to calculate climatic water balance for the area of Poland according to map algebra methodology (Fig. 5.5). Transformations included whole layers – as the variables all cells of the raster were used, their values were changed due to previously sited formulas for potential evapotranspiration and CWB index.

33 5.2. Validation

The last step was a validation of the proposed methods. The basis of spatial analysis is the comparison of real data and the data modelled with the suggested interpolation method at reference points. Because of the very few reference points as well as estimating different interpolation methods (which are not the typical methods of spatial data interpolation) only the common and unsophisticated statistical estimators were used (formulas 13-15). Through the approach the direct estimation of the model at each reference point. Within the last step of models validation also geostatistical estimators for spatial analyses were implemented (formulas 16-18).

The first and basis measure of model adjustment was the value of correlation coefficient between the real data (counted from in-situ measurements) and modelled data (R). The closer to value '1' an estimator value the better the model (closer to real image).

Also some basic statistical estimators (Stanisz 2007) presenting the difference between modelled and real values were implemented:

RE (real error) with equation:

(13)

where:

y_i – real value y’i – modelled value

PE (percentage error) with equation:

(14)

where:

yi – real value y’_i – modelled value

AE (absolute error) with equation:

| | (15)

where:

y_i – real value y’i – modelled value

34 In regards to a limited number of reference points it became impossible to estimate the models with the most common validation methods used for interpolation (neither cross- validation nor the method of independent sample could be properly used). However, estimators suggested to be used while analysing the results of spatial analyses conducted with different interpolation methods (ESRI 2001) were implemented to assess so called real spatial average interpolation errors, including info the formulas model errors for points gained in the first validation step:

RMSE (root-mean square error) with equation:

√^{∑ ( )} (16)

where:

y_i – real value y’i – modelled value n – pair points number

MPE (mean percentage error) with equation:

∑ (17)

where:

yi – real value y’i – modelled value n – pair points number

MAPE (mean absolute percentage error) with equation:

∑ | | (18)

where:

y_i – real value y’i – modelled value n – pair points number

The abovementioned estimators’ values closer to zero evidence the model quality.

All the model adjustment measures were implemented in relation to values in reference points.

35

6. S

P

Spatial differentiation of climatic water balance in Poland within growing season (May – October) received in the study amounted between less than -200 mm in central part of the country and hundreds millimetres in high mountain regions of the Carpathians Mountains and the Sudety Mountains (Fig 6.1).

Fig. 6.1. Spatial differentiation of climatic water balance (CWB in mm) in Poland in growing season (May–October) due to different methods: SLR – simple regression: f (H), MLR 1 – multiple regression: f (, H), MLR 2 – multiple regression: f (, H, d), MAG – map algebra

Most of the area of Poland is characterized by moisture shortage, positive values are usual only for the highlands, the forelands and the mountains areas (Fig 6.1). At the same time spatial distribution of the analysed index varies seasonally.

36 The highest deficiency is observed in May and June when spatial variability is also the most significant. Maximum values of CWB occur in October however their spatial distribution does not vary so much – almost the whole area of Poland is characterized by similar climatic water balance conditions. It is worth mentioning that spring months (May, June) are less humid than autumn ones (September, October) when CWB index posits higher values.

Nevertheless detailed analyses of CWB index distribution shows regional and local differences due to spatialisation method implemented which was described below.

37

6.1.

Regression models received in this study, mainly because of the significant correlation (compare with Table 5.2), are determined by altitude. The predominant role of this predictor has been influenced the spatial differentiation of climatic water balance index in Poland. Regardless the regression method, the beltlike distribution of CWB index values is less or more noticeable. This regularity is marked mainly in spring and summer months, not so visible in autumn (Figs. 6.2 – 6.4). Seasonal differentiation displays also in distance from the sea cost influence intensity. The impact of mentioned variable modifying climate conditions by evapotranspiration reduction and therefore higher CWB values is being reported from July to October while in May and June has no important effect.

Having the basic knowledge about spatial differentiation of climate conditions in Poland, it can be clearly stated that regression models give the CWB differentiation visualizations a bit deformed especially in previously described seaside area. Higher amount of water vapour in the atmosphere as well as higher precipitation totals do not result in moisture shortage within the range presented in Figures 6.2 – 6.4.

38 Fig. 6.2. Spatial differentiation of climatic water balance (CWB in mm) in Poland within growing season (May /V/–October/X/): SLR – simple regression: f (H)

39 Fig. 6.3. Spatial differentiation of climatic water balance (CWB in mm) in Poland within growing season (May /V/–October/X/): MLR 1 – multiple regression: f (, H)

40 Fig. 6.4. Spatial differentiation of climatic water balance (CWB in mm) in Poland within growing season (May /V/–October/X/): MLR 2 – multiple regression: f (, H, d)

41

6.2.

According to use of map algebra (MAG) climatic water balance field shows a significant moisture deficit (the lowest CWB index values) in lowlands area of Central Poland.

Fig. 6.5. Spatial differentiation of climatic water balance (CWB in mm) in Poland within growing season (May /V/–October/X/): MAG – map algebra

42 The situation is characteristic for all the growing period months (Fig. 6.5). Values well-marked higher than by other models appear at the cost, what confirm that the distance from the shore variable has been implemented in the MAG model. Beltlike distribution of CWB prove the index dependence on geographical parameters, e.g. altitude and mainly latitude.

Fundamental evidence acknowledge and supporting the MAG approach can be the detailed analysis of CWB components spatial fields which are the basis of the model (Fig.

6.6). The MAG image is the effect of spatialisation of different elements as precipitation (Fig. 4.10) and evapotranspiration (Fig. 6.9). Whereas the latter is a result of integrating solar radiation (Fig. 6.7) and temperature maps (Fig. 6.8).

Fig. 6.6. Spatial differentiation of CWB components in Poland in growing season (May–

October); A – solar radiation totals [cal·cm^-2]; B – temperature [°C], C – potential evapotranspiration [mm], D – precipitation totals [mm]

43 Fig. 6.7. Spatial differentiation of solar radiation [cal·cm^-2] in Poland within growing season (May /V/–October/X/)

44 Fig. 6.8. Spatial differentiation of air temperature [°C] in Poland within growing season (May /V/–October/X/)

45 Fig. 6.9. Spatial differentiation of potential evapotranspiration [mm] in Poland within growing season (May /V/–October/X/)

46 Fig. 6.10. Spatial differentiation of precipitation totals [mm] in Poland within growing season (May /V/–October/X/)

47

6.3.

Results achieved within both research methods were validated using universal statistical error estimators. CWB values calculated for 16 stations considered in the study accepted as the reference data.

For all the stations CWB values were calculated for the period under consideration (Table 5.1). For each point one defined the deviations of the modelled values from the real ones (CWB_mod – CWBcalc) within the growing season (May – October) and for each month separately. Positive error values (that is calculated differences) indicate model overestimation, whereas negative ones show undervaluation of the values predicted.

Table 6.1. Selected model errors of climatic water balance values in growing season (May – October)

errors MAG SLR MLR_1 MLR_2

MAE [mm] 48.6 82.0 83.0 75.5

RMSE [mm] 86.1 102.5 95.3 88.3

MAPE (%) 24.4 57.7 68.6 58.1

R 0.988 0.940 0.937 0.950

MAE – mean absolute error, RMSE – root-mean-square error, MAPE – mean absolute percentage error, R – Pearson’s correlation coefficient

MAG – map algebra, SLR – simple linear regression: f (H), MLR_1 – multiple linear regression: f (, H), MLR_2 – multiple linear regression: f (, H, d)

Validation results for the growing period (May – October) well-mark that map algebra (MAG) method is characterized by the smallest real errors (in relation to reference data) (Fig. 6.11, Table 6.1). The highest, close to 1, correlation coefficient value confirms the best model adjustment, the absolute errors are also significantly smaller than in the other methods (Table 6.1). Model overestimates values for the region of north and central Poland giving the best prediction for north-eastern part of the country and the worst one for the central part (Fig. 6.11). For southern Poland model values of climatic water balance are lower than the calculated (real). The differences reach at an average tens millimetres however, in extreme cases model can give CWB index value different from the real one of about few hundred millimetres. It has been reported at Kasprowy Wierch and Śnieżka mountains. Regression models which do not take into consideration the distance from the sea shore as the predictor (MLR 1) and simple regression (SLR) turned to be the weakest methods. Correlation coefficient about 0.5 lower and the other estimators achieve values a

48 little bit higher is in those cases (Table 6.1). It is worth noticing that spatial differentiation of the methods applicability is quite visible. Simple linear regression (SLR) gives better results for the areas of the Baltic coast and in the Sudety mountains (at lowlands differences between methods cannot be clearly conformed) whereas MLR 1 brings well- marked effects in the Carpathians (Fig. 6.11).

The worst results – regardless the method – were achieved for the mountain regions. All the models predict values lower that the real ones in the Carpathians (Table 6.1). For the Sudety mountains big positive differences, were modelled only for Śnieżka, the other stations are characterized by the errors similar to the rest of the country. The same difficult as in the case of mountains was to spatialise CWB index in the coastal area.

Regression models significantly reduce the prediction and estimate values close to those recorded for the lake districts and the central part of the country.

Some attempts were done to systematize the validation results for particular months. Every time the best results are given by map algebra (correlation coefficient values are not lower than 0.9), the worst ones while simple linear regression was implemented (correlation coefficient reaches less than 0.5). The most complicated was spatialisation of CWB index for May and June. The average of models absolute errors approximate 20 mm whereas in autumn such characteristic reaches about 15 mm (13 mm in October, 15 mm in September). It should be stated that high error values concerns especially regression models: SLR and MLR 1 (Tables 6.2., 6.3, 6.4 and 6.5). They correlation coefficient values are low and at the same time high values for other error estimators. The different situation takes place in autumn (in before mentioned September and October). Although regression models do not give the best results, the obtained fields can be prosperously compared with map algebra maps (Tables 6.6 and 6.7).

49 Fig. 6.11. Model errors [mm] and absolute percentage errors (%) of CWB in growing season at reference points (CWBmod–CWBref); SLR – simple regression: f (H), MLR 1 – multiple regression: f (, H), MLR 2 – multiple regression: f (, H, d), MAG – map algebra