Floodnet contribution

CSHShydRology.Rproj

@@ -15,4 +15,4 @@ LaTeX: pdfLaTeX
 BuildType: Package
 PackageUseDevtools: Yes
 PackageInstallArgs: --no-multiarch --with-keep.source
-PackageRoxygenize: rd,collate,namespace,vignette
+PackageRoxygenize: rd,collate,namespace

DESCRIPTION

@@ -21,11 +21,21 @@ Imports:
     timeDate,
     stringr,
     jsonlite,
-    curl
+    curl,
+    gstat,
+    sp,
+    lmom,
+    lmomco,
+    lmomRFA,
+    numDeriv,
+    mnormt,
+    Matrix,
+    boot,
+    quantreg
 Suggests: 
     knitr,
     testthat, 
     rmarkdown
 VignetteBuilder: knitr
 LazyData: true
-RoxygenNote: 6.1.1
+RoxygenNote: 7.0.2

NAMESPACE

@@ -1,5 +1,132 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(AIC,amax)
+S3method(AIC,fpot)
+S3method(AIC,nsamax)
+S3method(AIC,nsamax.mle)
+S3method(AIC,nspot)
+S3method(DataWide,data.frame)
+S3method(DataWide,formula)
+S3method(DataWide,matrix)
+S3method(DistSeason,data.frame)
+S3method(DistSeason,formula)
+S3method(DistSeason,matrix)
+S3method(DistSeason,numeric)
+S3method(ExtractAmax,default)
+S3method(ExtractAmax,formula)
+S3method(ExtractMmax,default)
+S3method(ExtractMmax,formula)
+S3method(FitPot,data.frame)
+S3method(FitPot,formula)
+S3method(FitPot,matrix)
+S3method(FitPot,numeric)
+S3method(GeoDist,default)
+S3method(GeoDist,formula)
+S3method(GofTest,amax)
+S3method(GofTest,fpot)
+S3method(GofTest,nspot)
+S3method(PadPot,data.frame)
+S3method(PadPot,default)
+S3method(PadPot,formula)
+S3method(RegSim,data.frame)
+S3method(RegSim,matrix)
+S3method(RegSim,reglmom)
+S3method(SeasonStat,data.frame)
+S3method(SeasonStat,default)
+S3method(SeasonStat,formula)
+S3method(coef,amax)
+S3method(coef,fpot)
+S3method(coef,nsamax)
+S3method(coef,nsamax.mle)
+S3method(coef,nspot)
+S3method(coef,poolmargin)
+S3method(coef,poolmle)
+S3method(coef,reglmom)
+S3method(fitted,nsamax)
+S3method(fitted,nspot)
+S3method(head,roicv)
+S3method(logLik,amax)
+S3method(logLik,fpot)
+S3method(logLik,nsamax)
+S3method(logLik,nsamax.mle)
+S3method(logLik,nspot)
+S3method(plot,amax)
+S3method(plot,fpot)
+S3method(plot,isite)
+S3method(plot,nsamax)
+S3method(plot,nspot)
+S3method(plot,reglmom)
+S3method(plot,roicv)
+S3method(predict,amax)
+S3method(predict,fpot)
+S3method(predict,nsamax)
+S3method(predict,nsamax.mle)
+S3method(predict,nspot)
+S3method(predict,poolmargin)
+S3method(predict,poolmle)
+S3method(predict,reglmom)
+S3method(predict,roi)
+S3method(print,amax)
+S3method(print,fpot)
+S3method(print,goftest)
+S3method(print,isite)
+S3method(print,nsamax)
+S3method(print,nsamax.mle)
+S3method(print,nspot)
+S3method(print,poolmargin)
+S3method(print,poolmle)
+S3method(print,reglmom)
+S3method(print,roi)
+S3method(residuals,nsamax)
+S3method(residuals,nspot)
+S3method(residuals,roi)
+S3method(simulate,amax)
+S3method(simulate,nsamax)
+S3method(simulate,nsamax.mle)
+S3method(simulate,nspot)
+S3method(simulate,poolmargin)
+S3method(simulate,poolmle)
+S3method(summary,bootns)
+S3method(vcov,amax)
+S3method(vcov,fpot)
+S3method(vcov,nsamax.mle)
+S3method(which.clusters,formula)
+S3method(which.clusters,numeric)
+S3method(which.floodPeaks,data.frame)
+S3method(which.floodPeaks,formula)
+S3method(which.floodPeaks,numeric)
+export(BootNsAmax)
+export(BootNsAmaxMle)
+export(BootNsPot)
+export(CvRoi)
+export(DataWide)
+export(DistSeason)
+export(ExtractAmax)
+export(ExtractMmax)
+export(FindNearest)
+export(FindThresh)
+export(FitAmax)
+export(FitNsAmax)
+export(FitNsAmaxMle)
+export(FitNsPot)
+export(FitPoolMargin)
+export(FitPoolMle)
+export(FitPot)
+export(FitRegLmom)
+export(FitRoi)
+export(GeoDist)
+export(GofTest)
+export(Intersite)
+export(JulianPlot)
+export(LmomDiag)
+export(PadPot)
+export(PlotMrl)
+export(PlotThresh)
+export(PoolRemove)
+export(RegSim)
+export(SearchThresh)
+export(SearchThreshNs)
+export(SeasonStat)
 export(ch_binned_MannWhitney)
 export(ch_booth_plot)
 export(ch_fdcurve)
@@ -18,16 +145,60 @@ export(ch_read_ECDE_flows)
 export(ch_regime_plot)
 export(ch_wtr_yr)
 export(cut_block)
+export(dAmax)
 export(date_subset)
+export(dgev)
+export(dglo)
+export(dgno)
+export(dgpa)
 export(doys)
+export(dpe3)
+export(fAmax)
+export(fgam)
+export(fgev)
+export(fglo)
+export(fgno)
+export(fgpa1d)
+export(fgpa2d)
+export(fgpaLmom)
+export(fgpaMom)
+export(fgum)
+export(fnor)
+export(fpe3)
+export(lAmax)
+export(npAmax)
+export(pAmax)
+export(pgev)
+export(pglo)
+export(pgno)
+export(pgpa)
+export(ppe3)
+export(qAmax)
+export(qgev)
+export(qglo)
+export(qgno)
+export(qgpa)
+export(qpe3)
+export(rAmax)
+export(rgev)
+export(rglo)
+export(rgno)
+export(rgpa)
+export(rpe3)
+export(sitenames)
 export(slice)
 export(stack_EC)
 export(sub_set_Years)
+export(which.clusters)
+export(which.day)
+export(which.floodPeaks)
 import(curl)
 import(fields)
 import(grDevices)
 import(graphics)
 import(jsonlite)
+import(lmom)
+import(lmomRFA)
 import(lubridate)
 import(stats)
 import(stringr)
@@ -41,6 +212,17 @@ importFrom(graphics,par)
 importFrom(graphics,plot)
 importFrom(graphics,points)
 importFrom(graphics,polygon)
+importFrom(lmomco,dlmomco)
+importFrom(lmomco,lmom2par)
+importFrom(lmomco,lmoms)
+importFrom(lmomco,mle2par)
+importFrom(lmomco,plmomco)
+importFrom(lmomco,pwm)
+importFrom(lmomco,pwm2lmom)
+importFrom(lmomco,qlmomco)
+importFrom(lmomco,rlmomco)
+importFrom(lmomco,vec2lmom)
+importFrom(lmomco,vec2par)
 importFrom(plotrix,radial.grid)
 importFrom(plotrix,radial.plot)
 importFrom(timeDate,as.timeDate)

R/rfa_AIC_amax.R

@@ -0,0 +1,3 @@
+#' @export
+AIC.amax <- function(object, k = 2, ...)
+  as.numeric(k*length(object$para) - 2*object$llik)

R/rfa_AIC_fpot.R

@@ -0,0 +1,4 @@
+#' @export
+#' @rdname FitPot
+AIC.fpot <- function(object, k = 2, ...)
+  -2 * logLik(object) + 2 * k

R/rfa_AIC_nsamax.R

@@ -0,0 +1,18 @@
+#' AIC for nonstationary AMAX model
+#' 
+#' Return the AIC of the nonstationary AMAX model. 
+#' 
+#' @author Martin Durocher <mduroche@@uwaterloo.ca>
+#'
+#' @param object Output
+#' @param ... Other parameters
+#' @param k Penalty per parameter.
+#'
+#' @export
+#'
+#' @seealso \link{AIC}.
+#' 
+AIC.nsamax <- function(object, ..., k = 2){
+  np <- length(coef(object))-1
+  return(k * np - 2*logLik(object)) 
+}

R/rfa_AIC_nsamax_mle.R

@@ -0,0 +1,6 @@
+#' @export
+#' @rdname AIC.nsamax
+AIC.nsamax.mle <- function(object, ..., k = 2){
+  np <- length(coef(object))
+  return(k * np - 2*logLik(object)) 
+}

R/rfa_BootNsAmax.R

@@ -0,0 +1,103 @@
+#' Bootstrap sample of the nonstationary AMAX model
+#' 
+#' Return bootstrap samples of the parameters and flood quantiles for the 
+#' nonstationary AMAX model.
+#' 
+#' @author Martin Durocher <mduroche@@uwaterloo.ca>
+#'
+#' @param object Output form \link{FitNsAmax}.
+#' 
+#' @param p Exceeding probability associated to a quantile or reliability level.
+#'    
+#' @param newdata Covariates. If \code{NULL} the input data that served to fit 
+#'   the model is used.
+#'   
+#' @param reliability Logical. Should the flood quantile associated with 
+#'   reliability level be returned. In that case, \code{newdata} must represent
+#'   the period of interest. See Durocher et al. (2019) for more details.
+#' 
+#' @param nsim Number of simulations.
+#' 
+#' @param alpha Confidence level.
+#' 
+#' @param verbose Logical. Should a progress bar be displayed.
+#'
+#' @param ... Other parameters.
+#' 
+#' @return 
+#' 
+#' \item{para}{Parameter of the detrended distribution.}
+#' \item{beta}{Parameter of the trend.}
+#' \item{qua}{flood quantiles.}
+#' 
+#' @seealso \link{FitNsAmax}, \link{predict.nsamax}.
+#'
+#' @export
+#'
+BootNsAmax <- 
+  function(object,
+           p = c(0.5, 0.8, 0.9, 0.95, 0.98, 0.99),
+           newdata = NULL,
+           reliability = FALSE,
+           nsim = 500, 
+           alpha = 0.95, 
+           verbose = TRUE){
+
+  if(is.null(newdata))
+    newdata <- object$data 
+
+  paras <- matrix(0, nsim, length(object$para$para))
+  betas <- matrix(0, nsim, length(object$beta))
+
+  if(reliability){
+    quas <- matrix(0, nsim, length(p))
+  } else {
+    quas <- array(0, dim = c(nrow(newdata), length(p), nsim))
+  }
+
+  ## Create a formula and data.frame for fitting simulation
+  xjj <- model.matrix(attr(object$data,'term'), object$data)
+  xjj <- data.frame(0, xjj)
+  newdata0 <- model.matrix(attr(object$data,'term'), newdata)
+  newdata0 <- data.frame(0, newdata0)
+  form0 <- as.formula(paste0(colnames(xjj)[1],'~ . -1'))
+
+  if(verbose)
+    bar <- txtProgressBar()
+
+  for(jj in 1:nsim){
+
+    if(verbose)
+      setTxtProgressBar(bar, jj / nsim)
+
+    ## Simulate a bootstrap sample
+    xjj[,1] <- simulate(object)
+
+    ## Fit the model
+    fit <- FitNsAmax(form = form0, x = xjj,
+                     distr = object$para$type, type = object$type) 
+
+    paras[jj,] <- fit$para$para
+    betas[jj,] <- fit$beta
+
+    if(reliability){
+      quas[jj,] <- predict(fit, p, newdata0, TRUE)
+    } else{
+      quas[,,jj] <- predict(fit, p, newdata0, FALSE)
+    }
+  }    
+
+  colnames(paras) <- names(object$para$para)
+  colnames(betas) <- names(object$beta)
+  colnames(quas) <- as.character(round(p,3))
+
+  ans <- list(para = paras,
+              beta = betas,
+              qua = quas)
+
+  class(ans) <- 'bootns'
+
+  return(ans)
+
+}
+