On 24/09/14 20:16, Martin Maechler wrote:
<SNIP>
> 1) has your proposal ever been provided in R?
> I'd be happy to add it to the robustX
> (http://cran.ch.r-project.org/web/packages/robustX) or even
> robustbase (http://cran.ch.r-project.org/web/packages/robustbase)
package.
<SNIP>
I have coded up the algorithm from the Cameron and Turner paper. Dunno
if it gives exactly the same results as my (Splus?) code from lo these
many years ago (the code that is lost in the mists of time), but it
*seems* to work.
It is not designed to work with actual "streaming" data --- I
don't know
how to do that. It takes a complete data vector as input. Someone who
knows about streaming data should be able to adapt it pretty easily.
Said he, the proverbial optimist.
The function code and a help file are attached. These files have had
their names changed to end in ".txt" so that they will get through the
mailing list processor without being stripped. With a bit of luck.
If they *don't* get through, anyone who is interested should contact me
and I will send them to you "privately".
cheers,
Rolf
--
Rolf Turner
Technical Editor ANZJS
-------------- next part --------------
rlas <- function(y,b=0.2,mfn=function(n){0.1*n^(-0.25)},
nstart=30,scon=NULL) {
# Initialize:
y0 <- y[1:nstart]
alpha <- median(y0)
s <- if(is.null(scon)) mean(abs(y0-alpha)) else scon
mu <- mfn(nstart)/s
eps <- s/nstart^b
kount <- sum(abs(alpha-y0) < eps)
g <- kount/(eps*nstart)
ny <- length(y)
n <- nstart+1
locn <- numeric(ny)
locn[1:(nstart-1)] <- NA
locn[nstart] <- alpha
scale <- numeric(ny)
scale[1:(nstart-1)] <- NA
scale[nstart] <- s
# Calculate recursively:
while(n <= ny) {
s <- if(is.null(scon)) ((n-1)*s + abs(y[n]-alpha))/n else scon
mu <- mfn(n)/s
eye <- if(abs(alpha - y[n]) < s/n^b) 1 else 0
g <- (1-1/n)*g + n^(b-1)*eye/s
a <- max(mu,g)
alpha <- alpha + sign(y[n]-alpha)/(a*n)
locn[n] <- alpha
scale[n] <- s
n <- n+1
}
list(locn=locn,scale=scale)
}
-------------- next part --------------
\name{rlas}
\alias{rlas}
\title{
Recursive location and scale.
}
\description{
Calculate a recursive estimate of location, asymptotically
equivalent to the median, and a recursive estimate of scale
equal to the \bold{MEAN} absolute deviation.
}
\usage{
rlas(y, b = 0.2, mfn = function(n){0.1 * n^(-0.25)},
nstart = 30, scon = NULL)
}
\arguments{
\item{y}{
A numeric vector of i.i.d. data whose location and scale
parameters are to be estimated.
}
\item{b}{
A tuning parameter (default value equal to that used by
Holst, 1987).
}
\item{mfn}{
Function of the index of the data which must be positive and
and tend to 0 as the index tends to infinity. The default
function is that used by Holst, 1987.
}
\item{nstart}{
Starting values for the algorithm are formed from the first
\code{nstart} values of \code{y}.
}
\item{scon}{
A constant value for the scale parameter \code{s}. If this
is provided then the algorithm becomes equivalent to the
algorithm of Holst, 1987.
}
}
\value{
A list with entries
\item{locn}{The successive recursive estimates of location. The
first \code{nstart - 1} of these are \code{NA}.}
\item{scale}{The successive recursive estimates of scale. If
\code{scon} is specified these are all equal to \code{scon}.}
}
\references{
Cameron, Murray A. and Turner, T. Rolf (1993). Recursive location and
scale estimators. \emph{Commun. Statist. --- Theory Meth.} \bold{22}
(9) 2503--2515.
Holst, U. (1987). Recursive estimators of location.
\emph{Commun. Statist. --- Theory Meth.} \bold{16} (8) 2201--2226.
}
\author{
\email{r.turner at auckland.ac.nz}
\url{http://www.stat.auckland.ac.nz/~rolf}
}
\examples{
set.seed(42)
y <- rnorm(10000)
z1 <- rlas(y)
z2 <- rlas(y,scon=1)
z3 <- rlas(y,scon=100)
OP <- par(mfrow=c(3,1))
plot(z1$locn,type="l")
abline(h=median(y),col="red")
plot(z2$locn,type="l")
abline(h=median(y),col="red")
plot(z3$locn,type="l")
abline(h=median(y),col="red")
par(OP)
}
}
\keyword{ univar }
\keyword{ robust }