s-news
[Top] [All Lists]

Re: [query]

To: "Adham Samia A Y" <samia.adham@ic.ac.uk>, "s-news@lists.biostat.wustl.edu" <s-news@lists.biostat.wustl.edu>
Subject: Re: [query]
From: ROGER NGOUENET <r.f.ngouenet@usa.net>
Date: 18 Dec 00 10:23:09 MST
Hi,

The following version of acf (acf.n - with a quick fix - ) tested only 
on S-plus 5.1 (Unix version) will add a title desired. 

e.g acf.n(rnorm(20), main = "test") generates the plot with test as title.

"acf.n" <- 
function(x, lag.max, type = "correlation", plot = T, main = NULL)
{
        ## Handle new and old style time series differently
        if(is(x, "series")) {
                isNewStyle <- T
                xOrig <- x
                xunit <- NULL
                freq <- 1
                x <- seriesData(x)
        }
        else {
                isNewStyle <- F
                x <- as.rts(x)
                tspx <- tspar(x)
                freq <- frequency(x)
                xunit <- units(x)
        }
        if(is.factor(x) || (is.data.frame(x) && any(sapply(x, "is.factor"))))
                stop("Cannot calculate the acf of factors")
        if(is.cts(x))
                freq <- 1
        if(is.data.frame(x))
                x <- numerical.matrix(x)
        else if(!is.matrix(x))
                x <- as.matrix(x)
        n <- dim(x)[1]
        nser <- dim(x)[2]
        dnames <- dimnames(x)[[2]]
        correlation <- partial <- FALSE
        itype <- charmatch(type, c("covariance", "correlation", "partial"),
                nomatch = 0)
        switch(itype + 1,
                stop("desired type of ACF is unknown"),
                {
                        correlation <- FALSE
                        type <- "covariance"
                }
                ,
                {
                        correlation <- TRUE
                        type <- "correlation"
                }
                ,
                {
                        partial <- TRUE
                        type <- "partial"
                }
                )
        x <- scale(x, T, F)
        storage.mode(x) <- "single"
        ier <- .Fortran("tsmiss",
                NAOK = T,
                x,
                as.integer(n),
                as.integer(nser),
                integer(1),
                integer(1),
                as.integer(0),
                specialsok = T)
        if(ier[[6]] > 0)
                stop("There are missing values in midst of time series")
        lo <- ier[[4]]
        ngood <- ier[[5]]
        x <- x[lo:(lo + ngood - 1),  ]
        if(missing(lag.max))
                lag.max <- as.integer(10 * log10(ngood/nser))
        else lag.max <- as.integer(lag.max)
        # make lag.max an integer
        if(lag.max <= 0) stop("lag.max must be positive")
        if(lag.max > ngood) {
                lag.max <- ngood - 1
                warning("lag.max > series length : reduced to series length - 1"
                        )
        }
        r <- array(0, dim = c(lag.max + 1, nser, nser))
        lag <- r
        storage.mode(r) <- "single"
        mlag <- .Fortran("mlag",
                as.integer(ngood),
                as.integer(ngood),
                as.integer(nser),
                as.integer(lag.max),
                x,
                single(nser),
                r = r)
        r <- mlag$r
        if(!partial)
                acf <- (.Fortran("mlag1",
                        as.integer(ngood),
                        as.integer(nser),
                        as.integer(lag.max),
                        r = r,
                        as.logical(correlation)))[[4]]
        else {
                a <- array(0, dim = c(nser, nser, lag.max))
                storage.mode(a) <- "single"
                partialacf <- a
                den <- array(0, dim = c(nser, nser))
                storage.mode(den) <- "single"
                acf <- .Fortran("ar",
                        as.integer(ngood),
                        as.integer(nser),
                        as.integer(lag.max),
                        r,
                        integer(1),
                        den,
                        a,
                        single(lag.max + 1),
                        single(lag.max + 1),
                        a,
                        a,
                        den,
                        den,
                        den,
                        den,
                        partialacf = partialacf,
                        ier = as.integer(0))
                if(acf$ier != 0)
                        stop("cannot invert matrix")
                acf <- aperm(acf$partialacf, c(3, 1, 2))
        }
        lag[1:(lag.max + 1),  ,  ] <- (0:lag.max)/freq
        for(i in 1:nser) {
                if(nser == i)
                        break
                for(j in (i + 1):nser) {
                        lag[, i, j] <- (0:lag.max)/freq
                        lag[, j, i] <-  - (0:lag.max)/freq
                }
        }
        if(partial)
                lag <- array(lag[-1,  ,  ], dim = c(lag.max, nser, nser))
        series.name <- deparse(substitute(x))
        z <- list(acf = acf, lag = lag, n.used = ngood, type = type, series = 
                series.name, units = xunit)
        if(!is.null(dnames))
                z$snames <- dnames
        if(plot) {
                acf.plot(z, if(nser <= 5) c(nser, nser) else c(5, 5), conf.int
                         = correlation || partial, main = main)
        }
        z
}












"Adham, Samia A Y" <samia.adham@ic.ac.uk> wrote:
Hi,

I wonder if anyone could help me to get rid of the title

"Series x" which appears automatically in the plot of the

"acf(x)" in S-Plus (Unix Version).


Thanks in advance,



Samia Adham 

PhD. student

Department of Mathematics (Statistics Section)

Imperial College, London
---------------------------------------------------------------------
This message was distributed by s-news@lists.biostat.wustl.edu.  To
unsubscribe send e-mail to s-news-request@lists.biostat.wustl.edu with
the BODY of the message:  unsubscribe s-news


____________________________________________________________________
Get free email and a permanent address at http://www.netaddress.com/?N=1

<Prev in Thread] Current Thread [Next in Thread>
  • Re: [query], ROGER NGOUENET <=