R/irreg.fd.R

Defines functions plot.sparse.fd irreg.fd

Documented in irreg.fd plot.sparse.fd

#' Sample Irregular Functional Data
#' @param mu function, scalar or a vector defining the mean function; default value: \code{0}.
#' @param X centered stochastic process defined by a function of the form 
#'          \code{X(tObs,n)} and returning \code{n*length(tObs)} matrix, where each row represents observations from a trajectory. Default value: \code{wiener.process()}.
#' @param n sample size; default value: \code{100}.
#' @param m a vector of sampling rate or scalar of average sampling rate or a function of the form \code{f(n)} generating \code{n} positive integers; default value: \code{5}.
#' @param sig standard deviation of measurement errors; if \code{NULL} then determined by \code{snr}.
#' @param snr  signal to noise ratio to determine \code{sig}; default value: \code{5}.
#' @param domain the domain; default value: \code{c(0,1)}.
#' @param delta the proportion of the domain to be observed for each trajectory; default value: \code{1}.
#' @details The number of observation for each trajectory is randomly generated by \code{rpois(m)+1}. For each trajectory, the reference time \code{Oi} is uniformly sampled from the interval \code{[domain[1]+delta*L/2,domain[2]-delta*L/2]}, where \code{L} is the length of \code{domain}, and the design points for the trajectory is uniformly sampled from the interval \code{[Oi-delta*L/2,Oi+delta*L/2]}.
#' @return a list with the following members 
#'      \describe{
#'          \item{\code{t}}{list of design points sorted in increasing order for each trajectory.}
#'          \item{\code{y}}{list of vectors of observations for each trajectory.}
#'      }
#'  and with attributes \code{sig}, \code{snr}, \code{domain}, \code{delta} and
#'      \describe{
#'          \item{y0}{\code{n*m} matrix of observations without measurement errors.}
#'      }    
#'      
#' @references 
#' \insertRef{Lin2020}{synfd}
#' 
#' @examples
#' # Gaussian trajectories with constant mean function 1
#' Y <- irreg.fd(mu=1, X=gaussian.process(), n=10, m=5)
#'
#' # trajectories froma a process defined via K-L representation
#' Y <- irreg.fd(mu=cos, X=kl.process(eigen.functions='FOURIER',distribution='LAPLACE'),n=10, m=5)
#'
#' # trajectories with specified individual sampling rate
#' Y <- irreg.fd(mu=1, X=gaussian.process(cov=matern), n=10, m=rpois(10,3)+2)
#' @export

irreg.fd <- function(mu=0, X=wiener.process(), n=100, m=5,
                      sig=NULL, snr=5, domain=c(0,1),
                      delta = 1)
{
    Lt <- list()
    Ly <- list()
    
    if(is.vector(m))
    {
        if(length(m) > 1)
            mi <- m
        else
            mi <- 1 + rpois(n,m)
    }
    else if(is.function(m))
    {
        mi <- m(n)
    }
    
    L <- domain[2]-domain[1]
    O <- runif(n,min=domain[1],max=domain[2]-delta*L)
    Lt <- sapply(1:n,function(i){
        s <- O[i]
        sort(runif(mi[i], min=s, max=s+delta*L))
    })
    
    
    # now construct Ly based on Lt
    Ly0 <- lapply(Lt, function(tobs){
        
        if (is.function(mu)) mui <- mu(tobs)
        else if(length(mu)==1) mui <- rep(mu,length(tobs))
        else stop('mu must be a scalar or a function.')
        y0 <- mui + X(tobs,1)
        return(y0)
    })
    
    if(!is.null(snr))
    {
        if(!is.infinite(snr))
        {
            # roughly estimate the expectation of the L2 norm of X
            EX2 <- mean(apply(X(seq(domain[1],domain[2],length.out=100),100)^2,
                              1,sum)*L/100
            )
            sig <- sqrt(EX2/snr)
        }
        else sig <- 0
    }
    
    if(sig != 0)
    {
        Ly <- lapply(Ly0,function(y){
            y + rnorm(n=length(y),sd=sig)
        })
    }
    else Ly <- Ly0
    
    R <- list(t=Lt,y=Ly)
    attr(R,'sig') <- sig
    attr(R,'snr') <- snr
    attr(R,'y0') <- Ly0
    attr(R,'domain') <- domain
    attr(R,'delta') <- delta
    attr(R,'class') <- 'sparse.fd'
    return(R)
}

#' Plot Irregular Functional Data
#' @param x the object generated by \code{\link{irreg.fd}}.
#' @param ... other parameters passed to \code{plot} and \code{lines}.
#' @importFrom graphics lines plot
#' @return a plot of the dataset \code{x}.
#' @export
plot.sparse.fd <- function(x,...)
{
    plot(x$t[[1]],x$y[[1]],...)
    if(length(x$t) > 1)
    {
        for(i in 2:length(x$t))
        {
            lines(x$t[[1]],x$y[[1]],...)
        }
    }
}

Try the synfd package in your browser

Any scripts or data that you put into this service are public.

synfd documentation built on July 1, 2020, 6:04 p.m.