Monthly time-series of water flow at Big Tree water gauge.

Description

Average monthly natural water flow (cubic feet per second) at the Big Tree gauge of the San Lorenzo river in Santa Cruz County, CA from 1937 through 2014.

Usage

BTflow

Format

A time series of length 936.

Source

https://waterdata.usgs.gov/nwis/

References

U.S. Geological Survey (2016). National Water Information System data available on the World Wide Web (USGS Water Data for the Nation). https://waterdata.usgs.gov/nwis/.

Daily time-series of ELI anomalies.

Description

ELI anomalies on the daily time-scale from January 1, 1979 to December 31, 2019 with all February 29ths omitted.

Usage

ELIanoms

Format

A time series of length 14965.

Source

https://portal.nersc.gov/archive/home/projects/cascade/www/ELI

References

Patricola, C.M., O’Brien, J.P., Risser, M.D. et al. Maximizing ENSO as a source of western US hydroclimate predictability. Clim Dyn 54, 351–372 (2020). doi: 10.1007/s00382-019-05004-8

Combines state space blocks of an exDQLM

Description

The function combines two models into a single state space model for an exDQLM.

Usage

combineMods(m1, m2)

Arguments

Value

List containing the new combined state space model components.

Examples

trend.comp = polytrendMod(2,rep(0,2),10*diag(2))
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))
model = combineMods(trend.comp,seas.comp)
# using dlm package
library(dlm)
model = combineMods(dlmModPoly(order=2,C0=10*diag(2)),dlmModTrig(365,2,C0=10*diag(4)))

Plot a component of an exDQLM

Description

The function plots the dynamic MAP estimates and 95% credible intervals (CrIs) of a specified component of an exDQLM. Alternatively, if just.theta=TRUE the MAP estimates and 95% credible intervals (CrIs) of a single element of the dynamic state vector are plotted.

Usage

compPlot(
  y,
  m1,
  index,
  add = FALSE,
  col = "purple",
  just.theta = FALSE,
  cr.percent = 0.95
)

Arguments

Value

A list of the following is returned:

• map.comp - MAP estimate of the dynamic component (or element of the state vector).

• lb.comp - Lower bound of the 95% CrIs of the dynamic component (or element of the state vector).

• ub.comp - Upper bound of the 95% CrIs of the dynamic component (or element of the state vector).

Examples

y = scIVTmag[1:365]
trend.comp = polytrendMod(2,rep(0,2),10*diag(2))
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))
model = combineMods(trend.comp,seas.comp)
M0 = exdqlmISVB(y,p0=0.85,model,df=c(0.98,1),dim.df = c(2,6),
                   gam.init=-3.5,sig.init=15,tol=0.05)
# plot first harmonic component
compPlot(y,M0,index=c(3,4),col="blue")

Create state space model of exDQLM from DLM

Description

The function creates a state space model of an exDQLM from "dlm" object.

Usage

dlmMod(m)

Arguments

Value

List containing only the components of m needed for the exDQLM state space model.

Examples

library(dlm)
m = dlmModPoly(order=2,C0=10*diag(2)) + dlmModTrig(365,2,C0=10*diag(4))
model = dlmMod(m)

exDQLM Diagnostics

Description

The function computes the following for the model(s) provided: the posterior predictive loss criterion based off the check loss, the one-step-ahead distribution sequence and its KL divergence from normality. The function also plots the following: the qq-plot and ACF plot corresponding to the one-step-ahead distribution sequence, and a time series plot of the MAP standard forecast errors.

Usage

exdqlmChecks(
  y,
  m1,
  m2 = NULL,
  plot = TRUE,
  cols = c("grey", "grey"),
  ref = NULL
)

Arguments

Value

A list containing the following is returned:

• m1.uts - The one-step-ahead distribution sequence of m1.

• m1.KL - The KL divergence of m1.uts and a standard normal.

• m1.pplc - The posterior predictive loss criterion of m1 based off the check loss function.

• m1.qq - The ordered pairs of the qq-plot comparing m1.uts with a standard normal distribution.

• m1.acf - The autocorrelations of m1.uts by lag.

If m2 is provided, analogous results for m2 are also included in the list.

Examples

y = scIVTmag[1:100]
model = polytrendMod(1,mean(y),10)
M0 = exdqlmISVB(y,p0=0.85,model,df=c(0.95),dim.df = c(1),
                  gam.init=-3.5,sig.init=15)
check.out = exdqlmChecks(y,M0,plot=FALSE)
check.out$m1.KL
check.out$m1

k-step-ahead Forecast

Description

The function estimates and plots the k-step-ahead forecasted quantile distribution from the filtered quantile estimates.

Usage

exdqlmForecast(
  y,
  start.t,
  k,
  m1,
  fFF = NULL,
  fGG = NULL,
  plot = TRUE,
  add = FALSE,
  cols = c("purple", "magenta"),
  cr.percent = 0.95
)

Arguments

Value

A list containing the following is returned:

• fa - The forecasted state mean vectors.

• fR - The forecasted state covariance matrices.

• ff - The forecasted quantile mean estimates.

• fQ - The forecasted quantile variances.

Examples

y = scIVTmag[1:100]
model = polytrendMod(1,quantile(y,0.85),10)
M0 = exdqlmISVB(y,p0=0.85,model,df=c(0.98),dim.df = c(1),
                   gam.init=-3.5,sig.init=15)
exdqlmForecast(y,start.t=90,k=10,M0)

exDQLM - ISVB algorithm

Description

The function applies an Importance Sampling Variational Bayes (ISVB) algorithm to estimate the posterior of an exDQLM.

Usage

exdqlmISVB(
  y,
  p0,
  model,
  df,
  dim.df,
  fix.gamma = FALSE,
  gam.init = NA,
  fix.sigma = TRUE,
  sig.init = NA,
  dqlm.ind = FALSE,
  exps0,
  tol = 0.1,
  n.IS = 500,
  n.samp = 200,
  PriorSigma = NULL,
  PriorGamma = NULL,
  verbose = TRUE
)

Arguments

Value

A list of the following is returned:

• run.time - Algorithm run time in seconds.

• iter - Number of iterations until convergence was reached.

• dqlm.ind - Logical value indicating whether gamma was fixed at 0, reducing the exDQLM to the special case of the DQLM.

• model - List of the state-space model including GG, FF, prior parameters m0 and C0.

• p0 - The quantile which was estimated.

• df - Discount factors used for each block.

• dim.df - Dimension used for each block of discount factors.

• sig.init - Initial value for sigma, or value at which sigma was fixed if fix.sigma=TRUE.

• seq.sigma - Sequence of sigma estimated by the algorithm until convergence.

• samp.theta - Posterior sample of the state vector variational distribution.

• samp.post.pred - Sample of the posterior predictive distributions.

• map.standard.forecast.errors - MAP standardized one-step-ahead forecast errors.

• samp.sigma - Posterior sample of scale parameter sigma variational distribution.

• samp.vts - Posterior sample of latent parameters, v_t, variational distributions.

• theta.out - List containing the variational distribution of the state vector including filtered distribution parameters (fm and fC) and smoothed distribution parameters (sm and sC).

• vts.out - List containing the variational distributions of latent parameters v_t.

If dqlm.ind=FALSE, the list also contains:

• gam.init - Initial value for gamma, or value at which gamma was fixed if fix.gamma=TRUE.

• seq.gamma - Sequence of gamma estimated by the algorithm until convergence.

• samp.gamma - Posterior sample of skewness parameter gamma variational distribution.

• samp.sts - Posterior sample of latent parameters, s_t, variational distributions.

• gammasig.out - List containing the IS estimate of the variational distribution of sigma and gamma.

• sts.out - List containing the variational distributions of latent parameters s_t.

Or if dqlm.ind=TRUE, the list also contains:

• sig.out - List containing the IS estimate of the variational distribution of sigma.

Examples

y = scIVTmag[1:1095]
trend.comp = polytrendMod(1,mean(y),10)
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))
model = combineMods(trend.comp,seas.comp)
M0 = exdqlmISVB(y,p0=0.85,model,df=c(1,1),dim.df = c(1,6),
                 gam.init=-3.5,sig.init=15,tol=0.05)

exDQLM - MCMC algorithm

Description

The function applies a Markov chain Monte Carlo (MCMC) algorithm to sample the posterior of an exDQLM.

Usage

exdqlmMCMC(
  y,
  p0,
  model,
  df,
  dim.df,
  fix.gamma = FALSE,
  gam.init = NA,
  fix.sigma = FALSE,
  sig.init = NA,
  dqlm.ind = FALSE,
  Sig.mh,
  joint.sample = FALSE,
  n.burn = 2000,
  n.mcmc = 1500,
  init.from.isvb = TRUE,
  PriorSigma = NULL,
  PriorGamma = NULL,
  verbose = TRUE
)

Arguments

Value

A list of the following is returned:

• run.time - Algorithm run time in seconds.

• model - List of the state-space model including GG, FF, prior parameters m0 and C0.

• p0 - The quantile which was estimated.

• df - Discount factors used for each block.

• dim.df - Dimension used for each block of discount factors.

• samp.theta - Posterior sample of the state vector.

• samp.post.pred - Sample of the posterior predictive distributions.

• map.standard.forecast.errors - MAP standardized one-step-ahead forecast errors.

• samp.sigma - Posterior sample of scale parameter sigma.

• samp.vts - Posterior sample of latent parameters, v_t.

• theta.out - List containing the distributions of the state vector including filtered distribution parameters (fm and fC) and smoothed distribution parameters (sm and sC).

If dqlm.ind=FALSE, the list also contains the following:

• samp.gamma - Posterior sample of skewness parameter gamma.

• samp.sts - Posterior sample of latent parameters, s_t.

• init.log.sigma - Burned samples of log sigma from the random walk MH joint sampling of sigma and gamma.

• init.logit.gamma - Burned samples of logit gamma from the random walk MH joint sampling of sigma and gamma.

• accept.rate - Acceptance rate of the MH step.

• Sig.mh - Covariance matrix used in MH step to jointly sample sigma and gamma.

Examples

y = scIVTmag[1:100]
trend.comp = polytrendMod(1,mean(y),10)
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))
model = combineMods(trend.comp,seas.comp)
M2 = exdqlmMCMC(y,p0=0.85,model,df=c(1,1),dim.df = c(1,6),
                gam.init=-3.5,sig.init=15,
                n.burn=100,n.mcmc=150)

Plot exDQLM

Description

The function plots the MAP estimates and 95% credible intervals (CrIs) of the dynamic quantile of an exDQLM.

Usage

exdqlmPlot(y, m1, add = FALSE, col = "purple", cr.percent = 0.95)

Arguments

Value

A list of the following is returned:

• map.quant - MAP estimate of the dynamic quantile.

• lb.quant - Lower bound of the 95% CrIs of the dynamic quantile.

• ub.quant - Upper bound of the 95% CrIs of the dynamic quantile.

Examples

y = scIVTmag[1:100]
model = polytrendMod(1,quantile(y,0.85),10)
M0 = exdqlmISVB(y,p0=0.85,model,df=c(0.98),dim.df = c(1),
                   gam.init=-3.5,sig.init=15)
exdqlmPlot(y,M0,col="blue")

Monthly Niño 3.4 Index.

Description

Monthly Niño 3.4 sea surface temperature (SST) Index time series.

Usage

nino34

Format

A time series of length 936.

Source

https://psl.noaa.gov/gcos_wgsp/Timeseries/Nino34/

References

Rayner N. A., D. E. Parker, E. B. Horton, C. K. Folland, L. V. Alexander, D. P. Rowell, E. C. Kent, A. Kaplan. Global analyses of sea surface temperature, sea ice, and night marine air temperature since the late nineteenth century. J. Geophys. Res., 108 (D14), 4407 (2003). doi: 10.1029/2002JD002670

Create an n-th order polynomial exDQLM component

Description

The function creates an n-th order polynomial exDQLM component.

Usage

polytrendMod(order, m0, C0)

Arguments

Value

A list of the following:

• FF - Observational vector.

• GG - Evolution matrix.

• m0 - Prior mean of the state vector.

• C0 - Prior covariance of the state vector.

Examples

# create a second order polynomial component
trend.comp = polytrendMod(2,rep(0,2),10*diag(2))

Time series of daily average magnitude IVT in Santa Cruz, CA.

Description

ECMWF Re-Analysis 5 (ERA5) daily average magnitude IVT in Santa Cruz, CA (approximately 22 N, 122 W) from January 1, 1979 to December 31, 2019 with all February 29ths omitted.

Usage

scIVTmag

Format

A time series of length 14965.

Source

https://cds.climate.copernicus.eu

References

Hersbach, H, Bell, B, Berrisford, P, et al. The ERA5 global reanalysis. Q J R Meteorol Soc. 2020; 146: 1999– 2049. doi: 10.1002/qj.3803

Create Fourier representation of a periodic exDQLM component

Description

The function creates a Fourier form periodic component for given period and harmonics.

Usage

seasMod(p, h, m0, C0)

Arguments

Value

A list of the following:

• FF - Observational vector.

• GG - Evolution matrix.

• m0 - Prior mean of the state vector.

• C0 - Prior covariance of the state vector.

Examples

# create a seasonal component with first, second and fourth harmonics of a period of 365
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))

Transfer Function exDQLM - ISVB algorithm

Description

The function applies an Importance Sampling Variational Bayes (ISVB) algorithm to estimate the posterior of an exDQLM with exponential decay transfer function component.

Usage

transfn_exdqlmISVB(
  y,
  p0,
  model,
  X,
  df,
  dim.df,
  lam,
  tf.df,
  fix.gamma = FALSE,
  gam.init = NA,
  fix.sigma = TRUE,
  sig.init = NA,
  dqlm.ind = FALSE,
  exps0,
  tol = 0.1,
  n.IS = 500,
  n.samp = 200,
  PriorSigma = NULL,
  PriorGamma = NULL,
  tf.m0 = rep(0, 2),
  tf.C0 = diag(1, 2),
  verbose = TRUE
)

Arguments

Value

A list of the following is returned:

• run.time - Algorithm run time in seconds.

• iter - Number of iterations until convergence was reached.

• dqlm.ind - Logical value indicating whether gamma was fixed at 0, reducing the exDQLM to the special case of the DQLM.

• model - List of the augmented state-space model including GG, FF, prior parameters m0 and C0.

• p0 - The quantile which was estimated.

• df - Discount factors used for each block, including transfer function component.

• dim.df - Dimension used for each block of discount factors, including transfer function component.

• lam - Transfer function rate parameter lambda.

• sig.init - Initial value for sigma, or value at which sigma was fixed if fix.sigma=TRUE.

• seq.sigma - Sequence of sigma estimated by the algorithm until convergence.

• samp.theta - Posterior sample of the state vector variational distribution.

• samp.post.pred - Sample of the posterior predictive distributions.

• map.standard.forecast.errors - MAP standardized one-step-ahead forecast errors.

• samp.sigma - Posterior sample of scale parameter sigma variational distribution.

• samp.vts - Posterior sample of latent parameters, v_t, variational distributions.

• theta.out - List containing the variational distribution of the state vector including filtered distribution parameters (fm and fC) and smoothed distribution parameters (sm and sC).

• vts.out - List containing the variational distributions of latent parameters v_t.

• median.kt - Median number of time steps until the effect of X_t is less than or equal to 1e-3.

If dqlm.ind=FALSE, the list also contains:

• gam.init - Initial value for gamma, or value at which gamma was fixed if fix.gamma=TRUE.

• seq.gamma - Sequence of gamma estimated by the algorithm until convergence.

• samp.gamma - Posterior sample of skewness parameter gamma variational distribution.

• samp.sts - Posterior sample of latent parameters, s_t, variational distributions.

• gammasig.out - List containing the IS estimate of the variational distribution of sigma and gamma.

• sts.out - List containing the variational distributions of latent parameters s_t.

Or if dqlm.ind=TRUE, the list also contains:

• sig.out - List containing the IS estimate of the variational distribution of sigma.

Examples

y = scIVTmag[1:1095]
X = ELIanoms[1:1095]
trend.comp = polytrendMod(1,mean(y),10)
seas.comp = seasMod(365,c(1,2,4),C0=10*diag(6))
model = combineMods(trend.comp,seas.comp)
M1 = transfn_exdqlmISVB(y,p0=0.85,model=model,
                          X,df=c(1,1),dim.df = c(1,6),
                          gam.init=-3.5,sig.init=15,
                          lam=0.38,tf.df=c(0.97,0.97))