[R] Re: L-moments

[HOSKING@watson.ibm.com]

> Dear R-users:
>
> do you know of an efficient algorithm for calculating sample L-moments
> (especially of higher order, 3 and up) under R (and/or Splus)?
> Any advice/help would be greatly appreciated.
>
> Janusz Kawczak.

Here is an R/Splus implementation of routine SAMLMU from my LMOMENTS
Fortran package (http://lib.stat.cmu.edu/general/lmoments).  It has
been lightly tested, but carries no guarantees of accuracy.

samlmr<-function(x,nmom=4){
  xok<-x!is.na(x)Ù
  n<-length(xok)
  if (nmom<=0) return(numeric(0))
  if (nmom<=2) rnames <- paste("l",1:nmom,sep="_")
  else rnames <- c("l_1","l_2",paste("t",3:nmom,sep="_"))
  lmom <- rep(NA,nmom)
  names(lmom) <- rnames
  if (n==0) return(lmom)
  xok<-sort(xok)
  nmom.actual<-min(nmom,n)
  lmom1Ù <- mean(xok)
  if (nmom.actual==1) return(lmom)
  temp<-seq(1-n,n-1,by=2)
  p1<-rep(1,n)
  p<-temp/(n-1)
  lmom2Ù<-mean(xok*p)
  if (nmom.actual==2) return(lmom)
  if (xok1Ù==xoknÙ) {
    warning("all data values equal")
    return(lmom)
  }
  for(j in 3:nmom.actual) {
    p2<-p1
    p1<-p
    p<-((2*j-3)*temp*p1-(j-2)*(n+j-2)*p2)/((j-1)*(n-j+1))
    lmomjÙ<-mean(xok*p)/lmom2Ù
  }
  return(lmom)
}

J. R. M. Hosking
hosking at watson.ibm.com
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._