| Title: | Stochastic Rainfall Generator Bartlett-Lewis Rectangular Pulse Model |
|---|---|
| Description: | Due to a limited availability of observed high-resolution precipitation records with adequate length, simulations with stochastic precipitation models are used to generate series for subsequent studies [e.g. Khaliq and Cunmae, 1996, <doi:10.1016/0022-1694(95)02894-3>, Vandenberghe et al., 2011, <doi:10.1029/2009WR008388>]. This package contains an R implementation of the original Bartlett-Lewis rectangular pulse model (BLRPM), developed by Rodriguez-Iturbe et al. (1987) <doi:10.1098/rspa.1987.0039>. It contains a function for simulating a precipitation time series based on storms and cells generated by the model with given or estimated model parameters. Additionally BLRPM parameters can be estimated from a given or simulated precipitation time series. The model simulations can be plotted in a three-layer plot including an overview of generated storms and cells by the model (which can also be plotted individually), a continuous step-function and a discrete precipitation time series at a chosen aggregation level. |
| Authors: | Christoph Ritschel |
| Maintainer: | Christoph Ritschel <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 1.0 |
| Built: | 2024-10-31 06:21:53 UTC |
| Source: | CRAN |
Beta.fun is a help function for OF
Beta.fun(a, b)Beta.fun(a, b)
a |
|
b |
|
beta returns value of Beta.fun for parameters a and b
Christoph Ritschel [email protected]
BL.acc accumulates the BLRPM stepfunction for a given accumulation time t.acc
at a given accumulation level acc.val. An offse can be defined. The unit is typically hours.
BL.acc(sfn, t.acc = 240, acc.val = 1, offset = 0)BL.acc(sfn, t.acc = 240, acc.val = 1, offset = 0)
sfn |
stepfunction of precipitation |
t.acc |
|
acc.val |
|
offset |
|
p.acc data.frame
Christoph Ritschel [email protected]
lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim acc.val <- 1 offset <- 0 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim) sfn <- BL.stepfun(simulation$cells) ts <- BL.acc(sfn,t.acc,acc.val,offset)lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim acc.val <- 1 offset <- 0 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim) sfn <- BL.stepfun(simulation$cells) ts <- BL.acc(sfn,t.acc,acc.val,offset)
Bartlett-Lewis Rectangular Pulse Model
BL.sim(lambda = 4/240, gamma = 1/10, beta = 0.3, eta = 2, mux = 4, t.sim = 240)BL.sim(lambda = 4/240, gamma = 1/10, beta = 0.3, eta = 2, mux = 4, t.sim = 240)
lambda |
value specifying the generation rate of storms [1/h] |
gamma |
value specifying the storm duration [1/h] |
beta |
value specifying the generation rate of cells [1/h] |
eta |
value specifying the cell duration [1/h] |
mux |
value specifying the cell intensity [mm/h] |
t.sim |
value specifying the simulation time [h] |
Model description (Rodriguez-Iturbe et al., 1987):
The model is a combination of 2 poisson processes and simulates storms and cells. During the given simulation time storms are generated in a poisson process with rate lambda. Those storms are given a exponetially distributed duration with parameter gamma. During its duration the storm generates in a second poisson process cells with rate beta. The first cell has to be instantaneous at the time of the storm arrival. The cell duration is exponentially distributed with parameter eta. For the whole lifetime each cell is given a constant intensity which is exponentially distributed with parameter 1/mux.
Aggregation:
The intensities of all cells alive at time t are summed up for total precipitation at time t.
Parameter estimation:
The model parameters (lambda,gamma,beta,eta,mux) can be estimated from simulated or
observed precipitation time series using the method of moments. Certain moments, e.g. mean, variance
can be calculated from the time series at different aggregation levels. These moments
can also be calculated theoretically from model parameters. Both sets of statistics can be
compared in an objective function, similar to a squared error estimator. By numerical optimization
the model parameters can be tuned to match the time series characteristics.
BL.sim generates model realisations of storms and cells by using given model
parameters lambda,gamma,beta,eta,mux for a given simulation time t.sim
BL.sim returns storms; data.frame of all storms containing information about occurence time, end time and number of cells
BL.sim returns cells; data.frame of all cells containing information about occurence time, end time, intensity and storm index
Christoph Ritschel [email protected]
lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim)lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim)
BL.stepfun calculates a continous stepfunction of precipitation from
the data.frame cells
BL.stepfun(cells)BL.stepfun(cells)
cells |
|
sfn returns stepfunction of precipitation
Christoph Ritschel [email protected]
lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim) stepfun <- BL.stepfun(simulation$cells)lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 simulation <- BL.sim(lambda,gamma,beta,eta,mux,t.sim) stepfun <- BL.stepfun(simulation$cells)
BLRPM.class defines a new class for objects of type BLRPM containing the information about storms,
cells, stepfunction and the precipitation time series.
BLRPM.classBLRPM.class
An object of class R6ClassGenerator of length 24.
Christoph Ritschel [email protected]
BLRPM.est estimates the five Bartlett-Lewis rectangular pulse model parameters lambda,gamma,beta,eta,mux
for a given time series data. At first the time series statistics at given accumulation levels acc.vals
are calculated. These statistics are given over to the parameter estimation algorithm together with
parameter starting values par. An objective function O.Fun can be specified, default is BLRPM.OF.
In addition the weights for different statistics and accumulation levels weights.mean, weights.var, weights.cov, weights.pz
can be specified. For the BLRPM objective function the user can select the measure of distance between observation
and model with OF: =1 quadratic, =2: quad extended, =3: absolute, =4: abs extended.
A scale parameter controls different cases in the objective function for differences in the scale of duration
parameters gamma and eta.
If a debugging is wished, debug can be set to TRUE and a log file is created in working directory.
Several optim parameters can be also defined. For specifics see ?optim.
BLRPM.est(RR,acc.vals,pars.in,O.Fun, weights.mean,weights.var,weights.cov,weights.pz,OF,debug, scale,method,lower,upper,use.log,maxit,ndeps,trace)BLRPM.est(RR,acc.vals,pars.in,O.Fun, weights.mean,weights.var,weights.cov,weights.pz,OF,debug, scale,method,lower,upper,use.log,maxit,ndeps,trace)
RR |
|
acc.vals |
|
pars.in |
|
O.Fun |
|
weights.mean |
|
weights.var |
|
weights.cov |
|
weights.pz |
|
OF |
|
debug |
set |
scale |
|
method |
|
lower |
|
upper |
|
use.log |
|
maxit |
|
ndeps |
|
trace |
|
$est returns vector of estimated parameters lambda,gamma,beta,eta,mux
$conv returns value of convergence of optimization, see optim for details
$mess returns character message about optimization if using "L-BFGS-B" method
$Z returns value of objective function for estimated parameters
Christoph Ritschel [email protected]
t.sim=240 lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 pars <- c(lambda,gamma,beta,eta,mux) sim <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim) est <- BLRPM.est(sim$RR,pars.in=pars,method="BFGS",use.log=TRUE)t.sim=240 lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 pars <- c(lambda,gamma,beta,eta,mux) sim <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim) est <- BLRPM.est(sim$RR,pars.in=pars,method="BFGS",use.log=TRUE)
BLRPM.OF is the objective function used for parameter estimation of the BLRPM parameters.
Given a set of BLRPM parameters par this function calculates a set of model statistics at
given accumulation time steps acc.vals. These model statistics are compared with given
time series statistics stats in the objective function. The user is able to define weights
for each statistic (has to be the same length as statistics input vector). Option for debugging is given.
A scale parameter defines a criterium for which different kinds of model statistics are calculated.
This criterium is mainly based on the timescale difference between storm duration parameter gamma and cell duration
parameter eta.
If use.log is true, the objective function needs logarithmic input parameters. The value
of OF defines the kind of objective function to be used: 1= quadratic 2= quadratic extended 3= absolute 4= absolute extended.
BLRPM.OF(par, stats, acc.vals = c(1, 3, 12, 24), weights = rep(1, length(stats)), debug = FALSE, scale = 1, use.log = TRUE, OF = 2)BLRPM.OF(par, stats, acc.vals = c(1, 3, 12, 24), weights = rep(1, length(stats)), debug = FALSE, scale = 1, use.log = TRUE, OF = 2)
par |
|
stats |
|
acc.vals |
|
weights |
|
debug |
|
scale |
|
use.log |
|
OF |
|
Z returns value of objective function for input parameters and input statistics
Christoph Ritschel [email protected]
BLRPM.sim is the main function for simulating precipitation with the Bartlett-Lewis rectangular pulse model.
It generates storms and cells using the given five BLRPM parameters lambda, gamma, beta, eta, mux for a given
simulation time t.sim. The function BLRPM.sim then accumulates a precipitation time series of length
t.akk (typically the same as t.sim) with an accumulation time step interval from the generated
storms and cells. An offset can be used to delay the precipitation time series for initialization reasons.
BLRPM.sim returns a list of different variables and data.frames: Storms, Cells, Stepfun, Precip, time.
BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc,interval,offset)BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc,interval,offset)
lambda |
|
gamma |
|
beta |
|
eta |
|
mux |
|
t.sim |
|
t.acc |
|
interval |
|
offset |
|
$storms returns data.frame containing information about storms: start, end, number of cells
$cells returns data.frame containing information about cells: start, end, intensity, storm index
$sfn returns stepfunction used to accumulate precipitation time series
$RR returns vector of accumulated precipitation with time step interval [mm/interval]
$time returns vector of time steps [interval]
Christoph Ritschel [email protected]
lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim interval <- 1 offset <- 0 simulation <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc=t.sim,interval,offset)lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim interval <- 1 offset <- 0 simulation <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc=t.sim,interval,offset)
Delta.fun is a help function for OF
Delta.fun(kappa, MStrich)Delta.fun(kappa, MStrich)
kappa |
|
MStrich |
|
Delta returns value of Delta.fun for kappa and MStrich
Christoph Ritschel [email protected]
class BLRPMplot.BLRPM plots an object of class BLRPM returned by the function BLRPM.sim with an option to plot
either only the storms and cells or to additionally plot the stepfunction and the precipitation time series
in a multiframe plot.
## S3 method for class 'BLRPM' plot(x, ..., OSC = FALSE, start.time = NULL, end.time = NULL, legend = TRUE, c.axis = 1.5, c.lab = 1.5, c.legend = 1.5)## S3 method for class 'BLRPM' plot(x, ..., OSC = FALSE, start.time = NULL, end.time = NULL, legend = TRUE, c.axis = 1.5, c.lab = 1.5, c.legend = 1.5)
x |
|
... |
Arguments to be passed to methods, such as graphical parameters (see |
OSC |
|
start.time |
|
end.time |
|
legend |
|
c.axis |
|
c.lab |
|
c.legend |
|
Christoph Ritschel [email protected]
lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim interval <- 1 offset <- 0 simulation <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc=t.sim,interval,offset) plot(simulation,OSC=FALSE) plot(simulation,OSC=TRUE,start.time=1,end.time=24)lambda <- 4/240 gamma <- 1/10 beta <- 0.3 eta <- 2 mux <- 4 t.sim <- 240 t.acc <- t.sim interval <- 1 offset <- 0 simulation <- BLRPM.sim(lambda,gamma,beta,eta,mux,t.sim,t.acc=t.sim,interval,offset) plot(simulation,OSC=FALSE) plot(simulation,OSC=TRUE,start.time=1,end.time=24)
TS.acc accumulates a given time series x at a given accumulation level acc.val. Minimum value
for acc.val is 2 [unit time]
TS.acc(x,acc.val)TS.acc(x,acc.val)
x |
|
acc.val |
|
x.acc TS.acc returns a vector of an accumulated time series
Christoph Ritschel [email protected]
x <- rgamma(1000,1) x.2 <- TS.acc(x,acc.val=2)x <- rgamma(1000,1) x.2 <- TS.acc(x,acc.val=2)
TS.stats calculates statistics of a given time series x at given accumulation
levels acc.vals. The calculated statistics are the mean of the first accumulation level,
the variance, auto-covariance lag-1 and the probability of zero rainfall of all given accumulation
levels of the time series. These statistics are needed for estimating the BLRPM parameters.
TS.stats(x,acc.vals)TS.stats(x,acc.vals)
x |
|
acc.vals |
|
stats TS.stats returns a vector of statistics calculated at given accumulation levels
Christoph Ritschel [email protected]
time.series <- rgamma(1000,shape=1) statistics <- TS.stats(time.series,acc.vals=c(1,3,12,24))time.series <- rgamma(1000,shape=1) statistics <- TS.stats(time.series,acc.vals=c(1,3,12,24))