VS-Lite in pure Python
Expected time to run through: 5 mins
In this tutorial, we demonstrate how to run VS-Lite in pure Python with PyVSL.
[1]:
%load_ext autoreload
%autoreload 2
import PyVSL
import pandas as pd
import numpy as np
In the below, we check the function first, and let’s pay attention to the units of the input T and P.
[2]:
PyVSL.VSL?
Signature:
PyVSL.VSL(
syear,
eyear,
phi,
T,
P,
T1=8,
T2=23,
M1=0.01,
M2=0.05,
Mmax=0.76,
Mmin=0.01,
alph=0.093,
m_th=4.886,
mu_th=5.8,
rootd=1000,
M0=0.2,
substep=0,
I_0=1,
I_f=12,
hydroclim='P',
)
Docstring:
Translated from VSLite_v2_3.m - Simulate tree ring width index given monthly climate inputs.
Basic Usage:
trw = VSLite_v2_3(syear,eyear,phi,T1,T2,M1,M2,T,P)
gives just simulated tree ring as ouput.
[trw,gT,gM,gE,M] = VSLite_v2_3(syear,eyear,phi,T1,T2,M1,M2,T,P))
also includes growth response to temperature, growth response to soil
moisture, scaled insolation index, and soil moisture estimate in outputs.
Basic Inputs:
syear = start year of simulation.
eyear = end year of simulation.
phi = latitude of site (in degrees N)
T1 = scalar temperature threshold below which temp. growth response is zero (in deg. C)
T2 = scalar temperature threshold above which temp. growth response is one (in deg. C)
M1 = scalar soil moisture threshold below which moist. growth response is zero (in v/v)
M2 = scalar soil moisture threshold above which moist. growth response is one (in v/v)
(Note that optimal growth response parameters T1, T2, M1, M2 may be estimated
using code estimate_vslite_params_v2_3.m also freely available at
the NOAA NCDC Paleoclimatology software library.)
T = (12 x Nyrs) matrix of ordered mean monthly temperatures (in degEes C)
P = (12 x Nyrs) matrix of ordered accumulated monthly precipitation (in mm)
Advanced Inputs (must be specified as property/value pairs):
'lbparams': Parameters of the Leaky Bucket model of soil moisture.
These may be specified in an 8 x 1 vector in the following
order (otherwise the default values are read in):
Mmax: scalar maximum soil moisture content (in v/v),
default value is 0.76
Mmin: scalar minimum soil moisture (in v/v), default
value is 0.01
alph: scalar runoff parameter 1 (in inverse months),
default value is 0.093
m_th: scalar runoff parameter 3 (unitless), default
value is 4.886
mu_th: scalar runoff parameter 2 (unitless), default
value is 5.80
rootd: scalar root/"bucket" depth (in mm), default
value is 1000
M0: initial value for previous month's soil moisture at
t = 1 (in v/v), default value is 0.2
substep: logical 1 or 0; perform monthly substepping in
leaky bucket (1) or not (0)? Default value is 0.
'intwindow': Integration window. Which months' growth responses should
be intregrated to compute the annual ring-width index?
Specified as a 2 x 1 vector of integer values. Both
elements are given in integer number of months since January
(July) 1st of the current year in the Northern (Southern)
hemisphere, and specify the beginning and end of the integration
window, respectively. Defaults is [1 ; 12] (eg. integrate
response to climate over the corresponding calendar year,
assuming location is in the northern hemisphere).
'hydroclim': Value is a single character either taking value ['P'] or ['M'].
If ['M'], then 9th input is interpreted as an estimate of
soil moisture content (in v/v) rather than as precipitation.
Model default is to read in precipitation and use the CPC's
Leaky Bucket model of hydrology to estimate soil moisture,
however if soil moisture observations or a more sophisticated
estimate of moisture accounting for snow-related processes
is available, then using these data directly are recommended
(and will also speed up code).
For more detailed documentation, see:
1) Tolwinski-Ward et al., An efficient forward model of the climate
controls on interannual variation in tree-ring width, Climate Dynamics (2011)
DOI: 10.1007/s00382-010-0945-5
2) Tolwinski-Ward et al., Erratum to: An efficient forward model of the climate
controls on interannual variation in tree-ring width, Climate Dynamics (2011)
DOI: 10.1007/s00382-011-1062-9
3) Tolwinski-Ward et al., Bayesian parameter estimation and
interpretation for an intermediate model of tree-ring width, Clim. Past
(2013), DOI: 10.5194/cp-9-1-2013
4) Documentation available with the model at http://www.ncdc.noaa.gov/paleo/softlib/softlib.html
File: ~/Github/PyVSL/PyVSL/core.py
Type: function
Load test data
First of all, we load the test data (can be downloaded with this link), which contains T and P at a grid point of the CCSM4 past1000 simulation. Note that the unit for T is degree C and that for P is accumulated monthly precipitation (mm/month).
[3]:
data_dict = pd.read_pickle('./data/test_T_P.pkl')
time = data_dict['time']
T = data_dict['T'] # unit: degC
P = data_dict['P'] # unit: mm/month
Run VS-Lite in pure Python
Now we are ready to run the wrapper to generate simulated TRW.
[4]:
%%time
T1, T2 = 1, 15
M1, M2 = 0.01, 0.05
res = PyVSL.VSL(
850, 1850, 45, # the starting year and ending year of the input T & P, along with the latitude
T, P,
T1=T1, T2=T2, M1=M1, M2=M2, # parameters of the thresholds
)
print(res.keys())
TRW = res['trw'] # this is the simulated TRW
print(TRW)
dict_keys(['trw', 'gT', 'gM', 'gE', 'Gr', 'M', 'potEv', 'width', 'width_mean', 'width_std'])
[ 0.03243917 1.62155893 2.00388961 ... -0.873099 0.07933166
-1.2066474 ]
CPU times: user 306 ms, sys: 2.23 ms, total: 308 ms
Wall time: 308 ms
[ ]: