Tutorial on Truncated Distribution for the Zero Flow

[ghost]

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")

[scworland]

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]

I had a tutorial issue written way back. The x2xlo Rd of lmomco has another example I think. I think I have an example too embedded in restore-2018. I am away from computer at the moment.

-wha (iPhone)

On Feb 18, 2018, at 2:01 PM, Scott Worland notifications@github.com wrote:

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.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/scworland-usgs/restore-2018/issues/12#issuecomment-366542989, or mute the thread https://github.com/notifications/unsubscribe-auth/AE5JQLPBngWzev_d3eXPnUEdiTC6PBFzks5tWIGOgaJpZM4RA94x .

[scworland]

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

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]

Quantiles above the pplo. There is a probability remapping, which is why I have the f2flo, flo2f, x2xlo functions in lmomco. I am hammered on a TxDOT proposal today, yes holiday!

-wha (iPhone)

On Feb 18, 2018, at 4:15 PM, Scott Worland notifications@github.com wrote:

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

[image: image] https://user-images.githubusercontent.com/11577981/36357646-fe890866-14ce-11e8-8d4a-eccad35c1da5.png

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?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/scworland-usgs/restore-2018/issues/12#issuecomment-366552889, or mute the thread https://github.com/notifications/unsubscribe-auth/AE5JQMJbPQ6wiWCqvltmQVEw6HLxUcykks5tWKEKgaJpZM4RA94x .

[ghost]

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 world is the 0th-percent probability in the fitted distribution world. (The domain on which the L-moments were originally computed, which here means the nonzero flows.)

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")

William H. Asquith, Research Hydrologist U.S. Geological Survey, Science Building MS-1053 Texas Tech University, Lubbock, Texas 79409 806-742-3129 work; 806-392-4148 cell

On Sun, Feb 18, 2018 at 4:15 PM, Scott Worland notifications@github.com wrote:

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

[image: image] https://user-images.githubusercontent.com/11577981/36357646-fe890866-14ce-11e8-8d4a-eccad35c1da5.png

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?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/scworland-usgs/restore-2018/issues/12#issuecomment-366552889, or mute the thread https://github.com/notifications/unsubscribe-auth/AE5JQMJbPQ6wiWCqvltmQVEw6HLxUcykks5tWKEKgaJpZM4RA94x .