Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Tutorial on Truncated Distribution for the Zero Flow #12

Closed
ghost opened this issue Dec 13, 2017 · 5 comments
Closed

Tutorial on Truncated Distribution for the Zero Flow #12

ghost opened this issue Dec 13, 2017 · 5 comments

Comments

@ghost
Copy link

ghost commented Dec 13, 2017

restore-2018/scripts/Asquith/zeroflow_dist_tutorial.txt

library(lmomco)
08167000 Guadalupe River Comfort, 1950s
These are the L-moments of NONZERO flow.

lmr <- lmomco::vec2lmom(c(114.5072, 80.29483, 0.6901283, 0.5524655))
kap <- lmomco::lmom2par(lmr, type="kap") # estimate the kappa distribution
pplo <- 0.078 # probability of zero flow for the entire decade
FF <- seq(0.001,.999,by=.001) # convenience nonexceedance probabilities

To begin, naively plot the fitted distribution to the L-moments

plot(qnorm(FF), lmomco::qlmomco(FF, kap), type="l", log="y", xaxt="n",
     xlab="", ylab="QUANTILE, CFS")

but this curve is wrong and over estimates the low end! So now we plot the distribution via a conditional probability offset or truncation. The key is in f2flo() and the transform (f - pp)/(1 - pp), which maps a [pplo, 1] interval into a [0,1] for purposes of computing the quantiles. For example, if the probability of zero flow is 0.078 as in this example, then the first fitted quantile at 0.078+eps in the real work is the 0th probability in the fitted distribution world. (The domain on which the L-moments were computed.)

lines(qnorm(  lmomco::f2f(  FF, pp=pplo)),
      qlmomco(lmomco::f2flo(FF, pp=pplo), kap), col=2, lwd=3)
add.lmomco.axis(las=2, tcl=0.5, side.type="NPP", npp.as.aep=TRUE) # x-axis
mtext("Probability of zeroflow is at 0.078, RED CURVE IS CORRECT, BLACK IS WRONG")

screen shot 2017-12-13 at 12 14 05 pm

@scworland
Copy link
Owner

scworland commented Feb 18, 2018

Will,

How would I modify this to handle pplo?

sw_lmom2q <- function(lmoms,FF){
  moms <- vec2lmom(lmoms, lscale=FALSE) 
  par <- paraep4(moms, snap.tau4=TRUE) 
  if(is.null(par)){
    Q <- NA
  }else{
    Q <- qlmomco(FF, par)
  }
  return(Q)
}

"lmoms" is currently a vector of L1,T2,T3, and T4. pplo could either be appended to the lmom vector or passed in separately. This is my guess,

sw_lmom2q <- function(lmoms,FF){
  moms <- vec2lmom(lmoms, lscale=FALSE) 
  par <- paraep4(moms, snap.tau4=TRUE) 
  if(is.null(par)){
    Q <- NA
  }else{
    Q <- qlmomco(f2flo(FF,pp=pplo),par) 
  }
  return(Q)
}

@ghost
Copy link
Author

ghost commented Feb 18, 2018 via email

@scworland
Copy link
Owner

scworland commented Feb 18, 2018

No problem, I think I have figured it out. One other question. Which quantiles does the truncated function return? For example,

image

the number of results vary based on the value of pplo. Because we are estimating 27 quantiles, is it the last x number of quantiles? Like number [2] above with 22 values, is that quantiles 5:27? If so, i would be safe to left pad everything with zeros to ensure that every value of FF has an associated quantile.

@ghost
Copy link
Author

ghost commented Feb 19, 2018 via email

@ghost
Copy link
Author

ghost commented Feb 20, 2018 via email

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant