# Posts Tagged ‘water stress’

The impacts of drought depend on time-scale. On short time-scales, drought means dry soil. On long time-scales, it means dry rivers and empty reservoirs. A region may simultaneously experience dry conditions on one time-scale and wet conditions on another e.g. wet soil but low streamflow or visa versa.

Standardized Precipitation Index (SPI) is a widely used measure of drought which can be defined for any time-scale of interest. For any location, SPI is normally distributed with zero mean and unit standard deviation. Index values > 2 indicate exceptionally wet conditions for that location, values < -2 indicate exceptionally dry conditions for that location, etc. Historical precipitation is the only input needed to compute SPI.

Australia experienced drought between 2002 and 2007. The image below shows SPI computed for a location in the drought-prone Murray-Darling basin of New South Wales. The time-series run from Jan 1948 to Jan 2010 and the index was calculated for time-scales from 1 to 12 months. Precipitation data is from NCEP Reanalysis [1] in a 1.875° × 1.875° grid cell centred at 30°S 145°E.

The drought of 2002 to 2007 shows up very clearly. It was preceeded by a wet period between 2005 and 2001. While 2009 showed an episode of severe drought at short time-scales, SPI at was normal/wet at longer time-scales during 2009. Agricultural yields recovered.

## Calculating SPI-M

Empirical rainfall probability distributions are far from normal (gaussian) and often approximate a shifted gamma distribution. The empirical cumulative probability distributions are used to transform the rainfall time-series into time-series of percentile probabilities. A normally distributed precipitation index is found by pretending that these percentile probabilities derive from a standard cumulative normal distribution and inverting to find the index values.

This is simple in R. If the vector data contains rainfall infall data, then:

fit.cdf <- ecdf(data) cdfs <- sapply(data,fit.cdf) SPI <- qnorm(cdfs) 

Tha rainfall data are M-month moving averages (current and previous months). A separate index is calculated for each calendar month to remove seasonality. The R code used to compute SPI values (based in NCEP Reanalysis or other data sets such as GCPC) is here.

[1] The NCEP/NCAR 40-year reanalysis project, Bull. Amer. Meteor. Soc., 77, 437-470, 1996

Noted Added 11 October 2011: I have uploaded a slightly improved SPI R script here. The function getPrecOnTimescale(precipitation,k) takes a vector of monthly precipitation values and returns a k-month average (i.e current month and prior k-1 months). getSPIfromPrec(precip.k) takes k-month precipitation values and returns the corresponding vector of SPI values.

West of Ireland oak woodland - October 11

Flux tower eddy covariance instruments measure CO2 flux between an ecosystem and the atmosphere. Observed CO2 fluxes (also called NEE, Net Ecosystem Exchange) can be correlated with locally measured solar intensity, temperature, humidity etc. With enough data, a detailed picture of the local response of the ecosystem to environmental variables to be built up. Flux towers are one of the important tools for quantitative understanding the sensitivity of vegetation to climatic variability.

CO2 flux data  are noisy, to say the least. If you want to understand ecosystem reponse, you need a statistical model.  This post gives some simple ideas and R code showing how such a model can be developed. The example used is Harvard Forest [1], Massachussetts, a temperate deciduous broadleaf forest for which a large amount of flux tower data (collected by Munger & Wofsy) is available. It could equally be applied to tropical forest, grassland, cropland etc at one of the 200+ fluxnet sites worldwide. The model falls short of the full scientific rigour given in the references, but it can easily be developed into a more complete tool.

The main interest here is ecosystem response to sunlight i.e. photosynthesis. The light dependent part of the NEE is GEP (Gross Ecosystem Production), while the light independent part is RE (ecosystem respiration).

### Pre-processing

As always, the first steps are about getting the data into the right form. We use half-hourly data for 1991-2007  for Harvard Forest  i.e Ameriflux Level 2, site code USHa-1. This raw data has not been subject to any gap-filling statistical treatment. Here are the pre-processing steps:

1. download the annual data files and create a dataframe harvard using the rbind() combining all of the data
2. replace missing data (e.g. -9999) with “NA”
3. add a new fornight index, indicating to which fortnight of the year each row of data belongs
4. select a subset of the Ameriflux variables of interest. e.g. column headings c(“YEAR”,”DTIME”,”FN”,”TA”,”NEE”,”TS1″,”VPD”,”PAR_in”)
5. drop rows of data when one or more variables are undefined
6. group the data according to fortnights index i.e. create a list of 26 elements where each element is a dataframe containing all data sharing the same fortnight index.
7. drop “dark” data i.e. rows of data where PAR < 30.

The environmental variables retained at step 4 are TA (atmospheric temperature), TS1  (soil temperature), VPD (vapour pressure deficit) and PAR (photosynthetically active radiation). PAR is solar radiation in the range 0.4-0.7 microns which is energetic enough to contribute to photosynthesis. It has units μ mols s-1 m-2. (1 μ mole = 6 x 1017.)

Here is what the R code looks like:  #download the annual data files and create a dataframe harvard loc=file.path("ftp://cdiac.ornl.gov/pub/ameriflux/data/Level2/Sites_ByID/US-Ha1/with_gaps/usmaharv_1991_L2.csv"); # Ameriflux US_Ha1 download.file(loc,"temp.dat"); ha_1991=read.csv("temp.dat", header=FALSE,skip=20); ......................... loc=file.path("ftp://cdiac.ornl.gov/pub/ameriflux/data/Level2/Sites_ByID/US-Ha1/with_gaps/usmaharv_2007_L2.csv"); # Ameriflux US_Ha1 download.file(loc,"temp.dat"); ha_2007=read.csv("temp.dat", header=FALSE,skip=20); harvard <- rbind(ha_1991,ha_1992,ha_....._2006,ha_2007); # merged dataframes #replace missing data with "NA" temp <- (harvard[,] == -9999 | harvard[,] == -6999); #replace -9999 or -6999 with NAs harvard[temp]=NA; #add a new fortnight index lims <- seq(from=1,to=365,by=14); # fortnightly data lims[27] <- 367; # bound adjust, last fortnight has extra days fortnight <- function(d) {sum(ifelse(d<lims,0,1))}  # convert from DOY to fortnights harvard <- data.frame( "FN" = sapply(harvardUST < 0.2|is.na(harvardTA); y <- is.na(harvardSTS1); u <- is.na(harvardSPAR_in); w <- !(x|y|z|u|v|u_cut); #boolean. w is true when all variables are defined and U* > 0.2 m/s harvardS <- data.frame("YEAR"=harvardFN[w],"TA"=harvardNEE[w], "TS1" = harvardVPD[w],"PAR" = harvardlatex \mathrm{NEE} = \mathrm{R_E} -{{\alpha \beta \mathrm{PAR}} \over {\beta + \alpha \mathrm{PAR}}}latex \alpha \rightarrow \alpha \epsilonlatex \epsilon = {\lambda^2 \over {\mathrm{\left(TA-25 \right)}^2 + \lambda^2} } {\mu^2 \over {\mathrm{VPD}^2 + \mu^2}}latex \lambdalatex \mulatex \mu = 2latex \alpha^{-1}\$ corresponds to the number of  photosynthetically active photons  arriving at the forest canopy for every CO2 molecule trapped by photosynthesis. This number is about 14.

 A common feature of future climate scenarios is competition between increased CO2 availability and increased water stress layed out on a global scale. Clearly flux tower data are good place to look for these effects. The R script for the model and plots can be found here. References [1] Factors controlling CO2 exchange on timescales from hourly to decadal at Harvard Forest, Urbanski, S., C. Barford, S. Wofsy, C. Kucharik, E. Pyle, J. Budney, K. McKain, D. Fitzjarrald, M. Czikowsky, J. W. Munger http://www.agu.org/pubs/crossref/2007.../2006JG000293.shtml [2] See for example, An evaluation of models for partitioning eddy covariance-measured net ecosystem exchange into photosynthesis and respiration P. Stoy et al, http://www.geos.ed.ac.uk/homes/pstoy/Stoy_06b.pdf [3] More detailed models are give in An empirical model simulating diurnal and seasonal CO2 flux for diverse vegetation types and climate conditions, M. Saito, S. Maksyutov, R. Hirata, and A. D. Richardson. http //www.biogeosciences.net/6/585/2009/bg-6-585-2009.html & references. 
 October 12, 2009 By joe Carbon, Climate, ecosystem model, Forests ecosystem model eddy covariance Terrestrial Carbon water stress No Comments 
 

 
 Categories Categories Select Category Agriculture Carbon Climate ecosystem model Finance and Economics Forests Land Cover Maps R Statistics Uncategorized Wind Energy /* <![CDATA[ */ (function() { var dropdown = document.getElementById( "cat" ); function onCatChange() { if ( dropdown.options[ dropdown.selectedIndex ].value > 0 ) { location.href = "http://joewheatley.net/?cat=" + dropdown.options[ dropdown.selectedIndex ].value; } } dropdown.onchange = onCatChange; })(); /* ]]> */ Archives Archives Select Month January 2018 September 2017 May 2017 April 2017 November 2016 September 2016 May 2016 October 2015 August 2015 June 2015 May 2015 April 2015 January 2015 September 2014 July 2014 June 2014 February 2014 October 2013 August 2013 June 2013 April 2013 January 2013 November 2012 August 2012 July 2012 April 2012 January 2012 November 2011 August 2011 July 2011 March 2011 February 2011 October 2010 August 2010 July 2010 May 2010 March 2010 February 2010 January 2010 December 2009 November 2009 October 2009 September 2009 July 2009 June 2009 May 2009 April 2009 March 2009 Recent Posts Trees, volcanoes and climate Braess Paradox in AC Power Networks On the speed of ships Sea temperatures from HMS Beagle 
 
 
 Meta Log in Entries RSS Comments RSS WordPress.org This website uses cookies to improve your experience. We'll assume you're ok with this, but you can opt-out if you wish.Accept Read MorePrivacy & Cookies Policy //<![CDATA[ jQuery(document).ready(function() { cli_show_cookiebar({ settings: '{"animate_speed_hide":"500","animate_speed_show":"500","background":"#fff","border":"#444","border_on":true,"button_1_button_colour":"#000","button_1_button_hover":"#000000","button_1_link_colour":"#fff","button_1_as_button":true,"button_2_button_colour":"#333","button_2_button_hover":"#292929","button_2_link_colour":"#444","button_2_as_button":false,"font_family":"inherit","header_fix":false,"notify_animate_hide":true,"notify_animate_show":false,"notify_div_id":"#cookie-law-info-bar","notify_position_horizontal":"right","notify_position_vertical":"bottom","scroll_close":false,"scroll_close_reload":false,"showagain_tab":true,"showagain_background":"#fff","showagain_border":"#000","showagain_div_id":"#cookie-law-info-again","showagain_x_position":"100px","text":"#000","show_once_yn":false,"show_once":"10000"}' }); }); //]]> var gaJsHost = (("https:" == document.location.protocol) ? "https://ssl." : "http://www."); document.write(unescape("%3Cscript src='" + gaJsHost + "google-analytics.com/ga.js' type='text/javascript'%3E%3C/script%3E")); try { var pageTracker = _gat._getTracker("UA-8154580-3"); // Cookied already: pageTracker._trackPageview(); } catch(err) {} //<![CDATA[ MapPluginInit( /* Default map width */ 500, /* Default map height */ 300, /* Use rel instad of title? */ false); //]]>