\documentclass[a4paper]{article} \usepackage{graphics} % Packages to allow inclusion of graphics \usepackage[pdftex]{graphicx} \usepackage{fancyhdr} \usepackage[figuresright]{rotating} \usepackage{natbib} % Margini \setlength{\textwidth} {170mm} \setlength{\textheight}{240mm} %Altezza testo 227 mm %\setlength{\topmargin} {0.1mm} \setlength{\evensidemargin}{-5mm} %Margini per l'opzione twoside \setlength{\oddsidemargin} {-5mm} \setlength{\topmargin}{-10mm} \SweaveOpts{keep.source=TRUE} %\VignetteIndexEntry{Figures 11 and 12 in Griffis and Stedinger (2007)} \title{Figures 11 and 12 in Griffis and Stedinger (2007)} \author{Alberto Viglione} \date{} \begin{document} \maketitle At page 489 Griffis and Stedinger (2007) write: \\ To allow use of the LP3 distribution in future L-moment frequency analyses, simple expressions were developed for $\tau_4$ as a function of $\tau_3$ of the form $$\tau_4 = a + b \tau_3 + c \tau_3^2 + d \tau_3^3$$ For given $\tau_3$, with the coefficients in Table 3, the approximations yield values of $\tau_4$ accurate to within 0.008 over the range $\tau_{3(min)} \leq \tau_3 \leq 0.9$, wherein Table 3 specifies the minimum value $\tau_{3(min)}$ for each value of $\gamma_x$ with $\sigma_x = 0$. The following function uses the above equation and Table 3 in Griffis and Stedinger (2007): <<>>= tau4_LP3 <- function (tau3_LP3) { # the approximations yield values of tau4 accurate to # within 0.008 over the range tau3min <= tau3 <= 0.9 gammax <- c(-1.4, -1, -0.5, 0, 0.5, 1, 1.4) table3abcd <- matrix(c(0.0602, -0.1673, 0.8010, 0.2897, 0.0908, -0.1267, 0.7636, 0.2562, 0.1166, -0.0439, 0.6247, 0.2939, 0.1220, 0.0238, 0.6677, 0.1677, 0.1152, 0.0639, 0.7486, 0.0645, 0.1037, 0.0438, 0.9327, -0.0951, 0.0776, 0.0762, 0.9771, -0.1394), ncol=4, nrow=7, byrow=TRUE, dimnames=list(paste("gammax=", gammax, sep=""), c("a","b","c","d"))) tau3s <- matrix(c(rep(1, length(tau3_LP3)), tau3_LP3, tau3_LP3^2, tau3_LP3^3), nrow=4, ncol=length(tau3_LP3), byrow=TRUE, dimnames=list(c("1","tau3","tau3^2","tau3^3"), paste("tau3=", tau3_LP3, sep=""))) tau3min <- c(-0.2308, -0.1643, -0.0740, 0, 0.0774, 0.1701, 0.2366) %*% t(rep(1, length(tau3_LP3))) tau3max <- rep(0.9, 7) %*% t(rep(1, length(tau3_LP3))) tau3s2 <- rep(1,7) %*% t(tau3_LP3) keeptau4 <- (tau3s2 >= tau3min*0.9999) & (tau3s2 <= tau3max*1.0001) tau4s <- table3abcd %*% tau3s tau4s[!keeptau4] <- NA return(tau4s) } @ Figure 11 is then given by: <>= gammax <- c(-1.4, -1, -0.5, 0, 0.5, 1, 1.4) tau3 <- seq(-0.2, 1, by=0.1) tau4 <- tau4_LP3(tau3) plot(c(-0.4, 1), c(0, 1), type="n", xlab="L-skewness, tau3", ylab="L-Kurtosis, tau4") grid() for (i in 1:7) { lines(tau3, tau4[i,], type="b", lty=i, pch=i) } legend("topleft", legend=c(gammax, "OLB"), title=expression(gamma[x]), pch=c(1:7,NA), lty=c(1:7,1), lwd=c(rep(1,7),2)) curve((5*x^2 - 1)/4, add=TRUE, lwd=2) @ \noindent {\bf Fig. 11.} L-moment ratio diagram for the LP3 distribution as a function of log space skew (OLB represents the overall theoretical lower bound of the $\tau_3 - \tau_4$ space). \vspace{10 mm} For Figure 12 load the library <<>>= library(nsRFA) @ and do the following: <>= t3b=-0.2; t3t=0.9; t4b=-0.1; t4t=0.8 plot(c(t3b+0.05, t3t-0.05), c(t4b, t4t), type="n", xlab=expression(tau[3]), ylab=expression(tau[4]), main="") grid() tipi <- c(1,2,1,4,5,6) spessori <- c(1,1.3,1.1,1.3,1.1,1.1) colori <- c(1,1,"darkgrey",1,1,1) tau3 <- seq(0, 0.9, by=0.1) tau4 <- tau4_LP3(tau3) tau4r <- apply(tau4, 2, range, na.rm=TRUE) polygon(x=c(tau3, rev(tau3)), y=c(tau4r[1,], rev(tau4r[2,])), density=20, col="darkgrey", border="darkgrey", angle=-45) GPA <- function(x) 0.20196*x + 0.95924*x^2 - 0.20096*x^3 + 0.04061*x^4 curve(GPA, t3b, t3t, add=TRUE, lty=tipi[5], lwd=spessori[5]) GEV <- function(x) { 0.10701 + 0.1109*x + 0.84838*x^2 - 0.06669*x^3 + 0.00567*x^4 - 0.04208*x^5 + 0.03763*x^6 } curve(GEV, t3b, t3t, add=TRUE, lty=tipi[4], lwd=spessori[4]) GLO <- function(x) 0.16667 + 0.83333*x^2 curve(GLO, t3b, t3t, add=TRUE, lty=tipi[6], lwd=spessori[6]) LN3 <- function(x) 0.12282 + 0.77518*x^2 + 0.12279*x^4 - 0.13638*x^6 + 0.11368*x^8 curve(LN3, t3b, t3t, add=TRUE, lty=tipi[1], lwd=spessori[1]) PE3 <- function(x) 0.1224 + 0.30115*x^2 + 0.95812*x^4 - 0.57488*x^6 + 0.19383*x^8 curve(PE3, t3b, t3t, add=TRUE, lty=tipi[2], lwd=spessori[2]) points(0, 0, pch=3, cex=1.2) points(0, 0.1226, pch=2, cex=1.2) points(1/3, 1/6, pch=5, cex=1.2) points(0.1699, 0.1504, pch=6, cex=1.2) points(0, 1/6, pch=4, cex=1.2) curve((5*x^2 - 1)/4, t3b, t3t, add=TRUE, lwd=2) legend("bottomright", c("EXP", "EV1", "LOG", "NOR", "UNIF"), pch=c(5, 6, 4, 2, 3), bty="n") legend("topleft", legend=c("LN3","P3","LP3","GEV","GP","GL","OLB"), lty=c(tipi,1), lwd=c(spessori,2), col=c(colori,1), bty="n") @ \noindent {\bf Fig. 12.} L-moment ratio diagram including GLO, GEV, LN, generalized Pareto GPA , P3, Gumbel, normal, and LP3 (light gray region represents LP3 distribution with $|\gamma_x| \leq 1.414$; {\it I have not plotted the dark gray region with restricted values of $\sigma_x$; I added some other distributions}. \begin{thebibliography}{} \bibitem[Griffis and Stedinger, 2007]{GriffisStedinger2007} Griffis, V.W., and Stedinger, J.R. (2007). \newblock Log-Pearson Type 3 distribution and its application in Flood Frequency Analysis. 1: Distribution characteristics. \newblock {\em Journal of Hydrologic Engineering}, 12(5):482--491, DOI:10.1061/(ASCE)1084-0699(2007)12:5(482). \end{thebibliography} \end{document}