The Nitrogen Risk Assessment Model for Scotland (version2; NIRAMS II)

Contents

  1. Model background and development

  2. Model stucture and conceptualisation

  3. Input data

  4. Installation and usage

  5. Model output

  6. Issues and future developments

1. Model background and development

1.1. Original NIRAMS model

The original Nitrogen Risk Assessment Model for Scotland (NIRAMS) was developed at the Macaulay Institute (now the James Hutton Institute) during 2002/3. Development work was led by Dr. Sarah Dunn and a description of the model is provided by Dunn et al. (2004a and 2004b). One of the main motivations was to create a simple, Scotland-specific nitrate leaching model, similar to the MAGPIE model for England and Wales developed by Lord and Anthony (2006). The model was designed with the EU Nitrates Directive in mind and key requirements included being able to apply it at national scale and in areas with sparse observational data.

Physical process representations were kept very simple in order to reduce the number of parameters as far as possible. NIRAMS is therefore highly conceptualised, with most of its computational complexity coming from the fact it is fully spatially distributed (i.e. grid based) rather than because any of the individual calculations are complex. This makes the model unusual, because historically most developers of water quality models have refined and increased the process complexity first, before expanding the spatial domain. For example, most fully distributed water quality models (e.g. MIKE SHE) are extremely complex and highly parameterised, and even semi-distributed models like INCA and HYPE have much more complicated process representation than NIRAMS.

Because it is spatially distributed NIRAMS is relatively slow to run, and it could also be criticised for being overly simplistic. On the other hand, the grid-based outputs are extremely flexible, because results from a single model run can be analysed in different ways for different spatial units. Furthermore, outside of a handful of academic study catchments, we generally lack the observational data necessary to robustly calibrate and evaluate highly parameterised models, so the simplicity of NIRAMS may be an advantage for many practical applications.

1.2. NIRAMS II

The original model was coded in Avenue for ArcView 3.x. Although Avenue is still a useful language for grid computations, it was officially discontinued in 2002 and subsequent releases of ArcGIS have favoured first Visual Basic and, most recently, Python. Starting in 2011, the NIRAMS Avenue code was ported to Python 2.7 and, at the same time, the model’s handling of agricultural data was redesigned by switching from parish-level census summaries to more detailed field- and farm-scale information based on the Scottish Integrated Administration and Control System (IACS). This adds significantly to the complexity of pre-processing the model's input data, but has the advantage of providing more spatial detail and allowing e.g. identification of diffuse pollution “hotspots”. Other changes have included refining the nitrate leaching representation by incorporating processes from the STREAM-N model (Dunn et al., 2013) and enhancing the model's performance and temporal resolution. The dependency on proprietary software (ESRI ArcGIS) has also been removed, by re-coding all the spatial analysis operations in pure Python. Nevertheless, the basic structure and conceptualisation of NIRAMS II remain broadly the same as for the original model.

2. Model stucture and conceptualisation

The model is spatially distributed with a 1 km by 1km grid resolution and is best considered in two parts – a water balance module and a nitrate leaching module (see image below). The water balance module incorporates information representing climate, soil properties and land use patterns in order to predict runoff, which is then passed to the nitrate leaching module. This combines the water balance results with detailed information representing agricultural activities to estimate the amount of nitrate leaching from the soil at each time step. The model considers three possible flow pathways for the leached nitrate: overland flow (OF), shallow sub-surface flow (SSF) and deeper groundwater flow (GW). The key outputs are $1 \; km^2$ resolution grids representing the amount of water and nitrate following each flow pathway at each time step. These grids can be aggregated both spatially and temporally to give estimates of the total runoff and amount of nitrate leached within particular areas and time periods of interest.

Despite its relatively simple conceptualisation, NIRAMS II includes a number of uncertain parameters which must be calibrated if the model is to generate plausible output. Some of the parameters associated with the water balance component (for example, the factors used to partition drainage between groundwater and sub-surface flow pathways) have been pre-calibrated for Scotland and are fixed in the version of the model released here. Other parameters can be adjusted by the user to control the model's output.

The current version of NIRAMS II includes the following user-specified parameters:

Parameter Description Comment
$T_{snow}$ Temperature at which precipitation begins to fall as snow ($^\circ C$) Plausible values range from a few degrees above to a few degress below zero.
$T_{melt}$ Temperature at which the snowpack begins to melt ($^\circ C$) Plausible values range from zero to a few degrees above zero.
Degree-day factor Parameter controlling the rate of melting (mm of water equivalent) for temperatures above $T_{melt}$ ($mm / ^\circ C/ day$) Plausible values range from just above zero to 10 or more, depending on regional climate and geography
Organic N factor Proportion of organic nitrogen that becomes immediately available for leaching (dimensionless) Plausible values range from about 0.15 to perhaps 0.70, depending on the dominant form organic N (e.g. slurry versus manure)
Mineralisation parameter Parameter controlling mineralisation rate (dimensionless) Plausible values for Scotland seem to be around 0.1 to 0.2
Denitrification parameter Parameter controlling denitrification rate (dimensionless) Plausible values for Scotland seem to be around 0.01 to 0.02
N leaching parameter Parameter controlling nitrogen leaching rate (dimensionless) Plausible values for Scotland seem to be around 1.0

3. Input data

NIRAMS II makes use of the following input datasets:

In order to make information from the agricultural census compatible with the other datasets and the model’s structure, some pre-processing of the information was required. The IACS data is collected annually by the Scottish Government and provides detailed information on the crops grown in the majority of agricultural fields across Scotland. Gaps in the datasets from 2001 to 2010 inclusive were first patched using information from the Land Cover Map 2007, and then combined with the results of the British Survey of Fertiliser Practice to obtain estimates for the application rate of inorganic nitrogen fertiliser in each field. These results were then aggregated to 1 km resolution using area-weighted averaging to produce gridded estimates of the amount of inorganic nitrogen applied to each grid cell in each year.

The JAC data are also collected annually and provide information on the number, age and type of livestock owned by each business receiving agricultural subsidies. The amount of organic nitrogen excreted annually by each animal class was taken from manure planning documentation issued to farmers within Scotland's Nitrate Vulnerable Zones (NVZs), and these figures were used to estimate the total amount of organic nitrogen produced each year by each business. This was then distributed spatially over appropriate land classes at the business scale, using a rule set designed to be broadly compatible with the application limits currently in force within the NVZs. This information was also aggregated to 1 km resolution to create a gridded time series of organic nitrogen application.

The annual estimates of the amount of organic and inorganic nitrogen applied are distributed temporally using a set of idealised time series (below) that define, for a variety of crop classes, the length of the growing season, the amount of nitrogen uptake and the timing of fertiliser application. The area under each time series sums to 1, so the amount of nitrogen applied to a given crop on any particular day can be estimated by multiplying the value in the series by the total annual nitrate load for that crop type.

Finally, all the gridded datasets were resampled to the same resolution and aligned to the same national scale bounding box, meaning they could be stored as simple compressed arrays within an HDF5 data store. This provides convenient and efficient data access and avoids the need to repeatedly resample or realign the datasets during model simulations.

4. Installation and usage

4.1. Prerequisites

To run the NIRAMS II model directly, the user must first obtain the input HDF5 file containing gridded datasets of meterological time series, soil properites and agricultural activities. Unfortunatey this file cannot be made freely available, as several of the datasets included in it were obtained under a license that precludes redistribution. If you are interested in using the model for academic research of relevance to the Scottish Government, it may be possible to make this file available - please contact the James Hutton Institute in Aberdeen to discuss this. Without the input file the model will not run, although it may be possible to generate your own input file or to adapt the code to suit your needs.

Having obtained the input file, the next step is to install Python 2.7. If you're on Windows, I highly recommend WinPython, which is a portable Python distribution that will not interfere with anything else on your system. You will also need to install GDAL. This can be fiddly, but with WinPython all you need to do is download the apropriate Python "wheel" for your Python distribution and then install it using the WinPython control panel (note that, as mentioned in the previous link, you may also need to download and install VCredist SP1).

4.2. Installing NIRAMS II

The model can be installed by visiting the repository home page and clicking the "Download ZIP" button towards the right-hand edge of the screen. Extract the .zip archive to any location on your system, then open the WinPython command prompt (not the standard Windows command prompt) and change directory to the "Code" folder (the one containing setup.py) by typing e.g.

cd D:\nirams_ii\Code

Then type

python setup.py install

to install the NIRAMS II package.

You should also make a working copy of the Excel file nirams_ii_input_template.xlsx, which is located in the top level of the downloaded directory. This file is used to set user-specified options and parameters to configure each model run. Note: it is important that you only edit the "Value" columns in this spreadsheet; modifying any other cells may lead to errors when the model is run.

Once the installer has completed, you can delete the .zip archive and the extracted folder from your system. If you need to download a fresh copy of the input template at any point, it can be accessed here.

4.3. Running the model

NIRAMS II can be run in two modes: single_run and param_combos. These names also refer to the worksheets in nirams_ii_input_template.xlsx, which allow the user to set options for each run type.

The simplest way to run NIRAMS II is in single_run mode. Here the user must open the single_run sheet of the Excel template and fill-in appropriate entries in the Value column for each row. The rows specify, for example, the location of the input HDF5 data file, the time period of interest, the output file to be created etc. There are also options to define the physical parameters for the model run. As described above, NIRAMS II has 7 physical parameters, which can be adjusted to control the model's behavior. In single_run mode, the user must specify a single value for each of these physical parameters.

param_combos mode allows the user to explore multiple different parameter sets. The code provided does not allow full auto-calibration or uncertainty estimation (see section 6 for further discussion of this), but it does make it easier to explore a range of different model parameterisations. To run the model in this mode, the user must fill-in the param_combos sheet of the Excel template. Most of the information required is the same as for a single run, but this time the user can enter comma-separated lists of values for each of the 7 physical parameters. The code then calculates all unique combinations of the specified parameters and runs the model once for each parameter set. For example, if the user entered single values for the first 5 physical parameters, but 0.01, 0.02 for the Denitrification parameter and 1.0, 1.1 for the N leaching parameter, the model would be run four times in total, with the first 5 parameters kept constant and the last two varied as follows:

Run ID Denitrification parameter N leaching parameter
1 0.01 1.0
2 0.01 1.1
3 0.02 1.0
4 0.02 1.1

In addition, in param_combos mode the user can specify the number of processors available for distributing the model runs. When this number is greater than 1, simple parallelisation of the code is achieved using Python's multiprocessing module.

Once the user has completed the appropriate worksheet of the Excel template, running the model is as simple as creating a Python script such as

if __name__ == '__main__':  

    from niramsii import nirams

    in_xls = r'D:\nirams_ii_input_template.xlsx' # Path to Excel template

    nirams.run_model(in_xls, 'single_run')       # Run in 'single_run' mode

Alternatively, the last line of the above script could be replaced by

    nirams.run_model(in_xls, 'param_combos')     # Run in 'param_combos' mode

Note: When running NIRAMS II in param_combos mode, it is important to place all your code within an

if __name__ == '__main__':

block. This is due to the way the multiprocessing module works on Windows - see e.g. here for a more detailed explanation. Placing if __name__ == '__main__' at the top of your script should not affect the behaviour of your code, but without it the multiprocessing component of NIRAMS II will generate recurring errors.

5. Model output

The model generates annual, 1 km resolution grids for the specified bounding box representing, for each grid cell, the amount of nitrate-N leached (kg/ha), the total drainage (mm) and the sub-surface drainage (mm). The sub-surface drainage is the sum of the water following the shallow sub-surface and groundwater flow pathways; the total drainage includes overland flow as well (see the diagram in section 2). These grids can be ued to estimate nitrate loadings and concentrations from the surface to any region (e.g. a groundwater body) within the specified bounding box.

The output grids are written to an HDF5 file and, optionally, as a series of GeoTiffs. The HDF5 output is more compact and portable, but the GeoTiffs are easier to work with in e.g. ArcGIS, which makes them useful for quickly exploring the output and making maps. Note that to display the GeoTiffs in ArcGIS you will first need to calculate statistics, otherwise they will appear as blank "NoData" grids.

The HDF5 files are easy and efficient to manipulate using code and they can also be explored using the ViTables GUI, which is included with WinPython. To access ViTables, simply open the WinPython command prompt (not the standard Windows command prompt) and type

vitables

The HDF5 files are designed to be self-describing, in the sense that all the parameter values and user options for each model run are stored within the output file itself. The screenshot below shows an output HDF5 from using the model in param_combos mode, with the key properties for each run listed under the "User attributes" tab.

Individual datasets store additional properties, such as the cell size, spatial reference etc.

To illustrate the kind of output generated by the model, the maps below show some of the results used to inform the 2013 Nitrates Directive review in Scotland. The left-hand image shows the estimated runoff grid, while the right-hand one shows nitrate-N concentrations averaged over Scotland's bedrock groundwater bodies.

The output grids can also be used to generate histograms of e.g. nitrate concentration for groundwater bodies of interest e.g.

For further details of the way the model has been used in the past, please see the CREW report on the 2013 Nitrates Directive review.

6. Issues and future developments

6.1. Process representation

NIRAMS II is designed to estimate nitrate loadings to the groundwater from the base of the soil profile, but it is not a groundwater model. During the 2013 Nitrates Directive review in Scotland, a reasonable relationship was established between the simulated loadings and observed groundwater nitrate concentrations, but in general there is no reason to expect such a relationship to exist. In Scotland, most groundwater movement is dominated by rapid fracture flow and the vast majority of aquifers are classified as low productivity. Groundwater nitrate concentrations are therefore reasonably responsive, such that in many locations the nitrate load from the surface averaged over a 5-year period is a reasonably good predictor of the level of sub-surface contamination. This relationship, however, does not hold for all groundwater bodies in Scotland: some, for example, are overlain by impermeable clay layers, which protect them from agricultural pollution at the surface; others are associated with iron- or sulphur-rich bedrock geology, which may be associated with high rates of groundwater denitrification, meaning that a groundwater body can safely assimilate larger nitrate loads. These issues are discussed in more detail in the CREW report, and planned improvements for NIRAMS II include accounting for such effects.

Other possibilities for extending the model include incorporating a more explicit representation of the groundwater, for example by adding a "groundwater box". This would lead to a more complex, highly parameterised model, but with improved capabilities for simulating long-term groundwater responses. As a starting point, the recently developed Water Balance Model (link not yet active; repo not yet publically available) would provide a good foundation for a new version of the model, but there are no plans to take this idea forward at present.

6.2. Calibration and evaluation

A wide variety of methods exist for auto-calibrating and evaluating complex environmental models. A summary of commonly employed techniques is provided here. Given the current structure of the NIRAMS II codebase, it would be relatively easy to link the model to a Markov chain Monte Carlo (MCMC) calibration framework, such as that provided by emcee. This would offer a much more comprehensive and robust means of assessing parameter-related uncertainty and sensitivity than is currently afforded by running the model in param_combos mode (see section 4.3).

One issue with using NIRAMS II within e.g. an MCMC framework is that, when run at national scale, the model is slow. Although the present code is by no means highly optimised, it seems likely that uncertainty estimation procedures based on stochastic sampling will always be slow when run with spatially distributed models such as NIRAMS II (because evaluating the likelihood function is computationally "expensive"). One possible solution is to build a statistical emulator of the main model, which can then be used to explore the state space more rapidly using e.g. MCMC. This approach is currently being investigated by researchers at Biomathematics and Statistics Scotland and may form part of the upcoming 2017 Nitrates Directive review.