[Robkalman-commits] r30 - in branches/robKalman_bs/pkg: . robKalman robKalman/R robKalman/chm robKalman/demo robKalman/inst robKalman/man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Tue Jun 9 20:21:49 CEST 2009


Author: ruckdeschel
Date: 2009-06-09 20:21:49 +0200 (Tue, 09 Jun 2009)
New Revision: 30

Added:
   branches/robKalman_bs/pkg/robKalman/
   branches/robKalman_bs/pkg/robKalman/DESCRIPTION
   branches/robKalman_bs/pkg/robKalman/NAMESPACE
   branches/robKalman_bs/pkg/robKalman/R/
   branches/robKalman_bs/pkg/robKalman/R/0AllClass.R
   branches/robKalman_bs/pkg/robKalman/R/ACMfilt.R
   branches/robKalman_bs/pkg/robKalman/R/ACMfilter.R
   branches/robKalman_bs/pkg/robKalman/R/AllGeneric.R
   branches/robKalman_bs/pkg/robKalman/R/AllGenerics.R
   branches/robKalman_bs/pkg/robKalman/R/AllInitialize.R
   branches/robKalman_bs/pkg/robKalman/R/AllPlot.R
   branches/robKalman_bs/pkg/robKalman/R/AllShow.R
   branches/robKalman_bs/pkg/robKalman/R/Psi.R
   branches/robKalman_bs/pkg/robKalman/R/Util.R
   branches/robKalman_bs/pkg/robKalman/R/arGM.R
   branches/robKalman_bs/pkg/robKalman/R/arGMinternal.R
   branches/robKalman_bs/pkg/robKalman/R/calibrateRLS.R
   branches/robKalman_bs/pkg/robKalman/R/classKalman.R
   branches/robKalman_bs/pkg/robKalman/R/mACMfilter.R
   branches/robKalman_bs/pkg/robKalman/R/mACMinternal.R
   branches/robKalman_bs/pkg/robKalman/R/rLSfilter.R
   branches/robKalman_bs/pkg/robKalman/R/recFilter.R
   branches/robKalman_bs/pkg/robKalman/R/simulateSScont.R
   branches/robKalman_bs/pkg/robKalman/chm/
   branches/robKalman_bs/pkg/robKalman/chm/00Index.html
   branches/robKalman_bs/pkg/robKalman/chm/0robKalman-package.html
   branches/robKalman_bs/pkg/robKalman/chm/ACMfilt.html
   branches/robKalman_bs/pkg/robKalman/chm/Logo.jpg
   branches/robKalman_bs/pkg/robKalman/chm/Rchm.css
   branches/robKalman_bs/pkg/robKalman/chm/arGM.html
   branches/robKalman_bs/pkg/robKalman/chm/calibrateRLS.html
   branches/robKalman_bs/pkg/robKalman/chm/internalACM.html
   branches/robKalman_bs/pkg/robKalman/chm/internalKalman.html
   branches/robKalman_bs/pkg/robKalman/chm/internalarGM.html
   branches/robKalman_bs/pkg/robKalman/chm/internalpsi.html
   branches/robKalman_bs/pkg/robKalman/chm/internalrLS.html
   branches/robKalman_bs/pkg/robKalman/chm/recFilter.html
   branches/robKalman_bs/pkg/robKalman/chm/robKalman.chm
   branches/robKalman_bs/pkg/robKalman/chm/robKalman.hhp
   branches/robKalman_bs/pkg/robKalman/chm/robKalman.toc
   branches/robKalman_bs/pkg/robKalman/chm/simulateSScont.html
   branches/robKalman_bs/pkg/robKalman/chm/util.html
   branches/robKalman_bs/pkg/robKalman/demo/
   branches/robKalman_bs/pkg/robKalman/demo/00Index
   branches/robKalman_bs/pkg/robKalman/demo/ACMdemo.R
   branches/robKalman_bs/pkg/robKalman/demo/rLSdemo.R
   branches/robKalman_bs/pkg/robKalman/inst/
   branches/robKalman_bs/pkg/robKalman/inst/NEWS
   branches/robKalman_bs/pkg/robKalman/man/
   branches/robKalman_bs/pkg/robKalman/man/0robKalman-package.Rd
   branches/robKalman_bs/pkg/robKalman/man/ACMfilt.Rd
   branches/robKalman_bs/pkg/robKalman/man/arGM.Rd
   branches/robKalman_bs/pkg/robKalman/man/calibrateRLS.Rd
   branches/robKalman_bs/pkg/robKalman/man/internalACM.Rd
   branches/robKalman_bs/pkg/robKalman/man/internalKalman.Rd
   branches/robKalman_bs/pkg/robKalman/man/internalarGM.Rd
   branches/robKalman_bs/pkg/robKalman/man/internalpsi.Rd
   branches/robKalman_bs/pkg/robKalman/man/internalrLS.Rd
   branches/robKalman_bs/pkg/robKalman/man/recFilter.Rd
   branches/robKalman_bs/pkg/robKalman/man/simulateSScont.Rd
   branches/robKalman_bs/pkg/robKalman/man/util.Rd
Log:
@ Bernhard: arg-- hatte Ebene pkg vergessen... (Schritt 2)

Added: branches/robKalman_bs/pkg/robKalman/DESCRIPTION
===================================================================
--- branches/robKalman_bs/pkg/robKalman/DESCRIPTION	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/DESCRIPTION	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,14 @@
+Package: robKalman
+Version: 0.3
+Date: 2009-03-18
+Title: Robust Kalman Filtering
+Description: Routines for Robust Kalman Filtering --- the ACM- and rLS-filter
+Author: Peter Ruckdeschel, Bernhard Spangl, Irina Ursachi, Cezar Chirila
+Maintainer: Peter Ruckdeschel <Peter.Ruckdeschel at itwm.fraunhofer.de>
+Depends: R(>= 2.3.0), methods, graphics, startupmsg, dse1, dse2, MASS, limma, robustbase, numDeriv 
+Imports: stats, MASS
+LazyLoad: yes
+License: LGPL-3
+URL: http://robkalman.r-forge.r-project.org/
+LastChangedDate: {$LastChangedDate: 2009-03-31 15:31:30 +0200 (Di, 31 Mrz 2009) $}
+LastChangedRevision: {$LastChangedRevision: 447 $}

Added: branches/robKalman_bs/pkg/robKalman/NAMESPACE
===================================================================
--- branches/robKalman_bs/pkg/robKalman/NAMESPACE	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/NAMESPACE	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,8 @@
+import("methods")
+import("stats")
+import("startupmsg")
+
+export("ACMfilt", "ACMfilter", "arGM", "Euclidnorm",  
+       "simulateState", "simulateObs", "rcvmvnorm", "Huberize",
+       "rLScalibrateB", "limitS", "rLSFilter", "KalmanFilter", 
+       "recursiveFilter")

Added: branches/robKalman_bs/pkg/robKalman/R/0AllClass.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/0AllClass.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/0AllClass.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,68 @@
+.onLoad <- function(lib, pkg){
+    require("methods", character = TRUE, quietly = TRUE)
+
+}
+
+
+.onAttach <- function(library, pkg)
+{
+buildStartupMessage(pkg="robKalman", library=library, packageHelp=TRUE #, 
+                    #MANUAL=""
+                    )
+  invisible()
+}
+
+#.onUnload <- function(libpath)
+#{
+#    library.dynam.unload("distrEx", libpath)
+#}
+#
+#
+## register zoo as "S4"-class
+setOldClass("zoo")
+#
+setClassUnion("Hyperparamtype", 
+               c("NULL","matrix", "array", "zoo")
+               )
+
+#
+#
+#
+#
+## positive definite, symmetric matrices with finite entries
+#setClass("ACMcontrol",representation())
+#setClass("rLScontrol",representation())
+#
+#
+#
+#
+# class SSM --- State space model
+setClass("SSM",
+          representation = representation(
+                                name = "character",   ## name of the ssm
+                                F = "Hyperparamtype", ## transition matrix/ces or NULL
+                                Z = "Hyperparamtype", ## observation matrix/ces or NULL
+                                Q = "Hyperparamtype", ## innovation covariance or NULL
+                                V = "Hyperparamtype", ## observation error covariance or NULL
+                                p = "numeric",  ## state dimension
+                                q = "numeric"), ## observation dimension
+          prototype = prototype(name = gettext("a state space"), 
+                                F = NULL,
+                                Z = NULL,
+                                Q = NULL,
+                                V = NULL,
+                                p = 1, 
+                                q = 1), 
+          contains = "VIRTUAL")
+
+# class TimeInvariantSSM 
+setClass("TimeInvariantSSM",
+          prototype = prototype(name = gettext("a time-invariant state space"), 
+                                F = 1,
+                                Z = 1,
+                                Q = 1,
+                                V = 1,
+                                p = 1, 
+                                q = 1), 
+          contains = "SSM")          
+          
\ No newline at end of file

Added: branches/robKalman_bs/pkg/robKalman/R/ACMfilt.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/ACMfilt.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/ACMfilt.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,78 @@
+ACMfilt <- function (x, gm, s0=0, 
+                     psi="Hampel", a=2.5, b=a, c=5.0, 
+                     flag="weights", lagsmo=TRUE)
+{
+###########################################
+##
+##  R-function: ACMfilt - approximate conditional-mean filtering (wrapper)
+##  author: Bernhard Spangl
+##  version: 1.1 (2007-08-13 and 2006-08-31)
+##  References: 
+##  [Mart79c] R.D. Martin, Approximate Conditional-mean Type Smoothers 
+##                and Interpolators (1979)
+##  [Mart81b] R.D. Martin, Robust Methods for Time Series (1981)
+##  [MarT82b] R.D. Martin & D.J. Thomson, Robust-resistent Spectrum 
+##                Estimation (1982)
+##
+###########################################
+
+##  Paramters:
+##  x ... univariate time series (vector)
+##  gm ... list as produced by function 'arGM' which includes components 
+##         'ar' containing the AR(p) coefficient estimates, 'sinnov' containing 
+##         innovation scale estiamtes from AR(p) fits of orders 1 through p;
+##         'Cx' containing an estimate of the p by p autocovariance matrix, 
+##         and 'mu', the estimated mean of 'x'. 
+##  s0 ... scale of nominal Gaussian component of the additive noise
+##  psi ... influence function to be used (default: "Hampel", 
+##          only Hampel's psi function available at the moment)
+##  a, b, c ... tuning constants for Hampel's psi-function
+##              (defaul: a=b=2.5, c=5.0)
+##  flag ... character, if "weights" (default), use psi(t)/t to calculate 
+##           the weights; if "deriv", use psi'(t)
+##  lagsmo ... logical, if TRUE (default) lag p-1 smoothing is performed; 
+##             if FALSE filtering from the top of ^X_t is performed
+
+##  Variable definitions:
+
+    N <- length(x)
+    phi <- gm$ar
+    p <- length(phi)
+    si <- gm$sinnov[p]
+    Cx <- gm$Cx
+    Phi <- cbind(rbind(phi[-p], diag(rep(1, (p-1)))), c(phi[p], rep(0, (p-1))))
+    Q <- matrix(0, p, p)
+##  Q <- diag(rep(0, p))
+    Q[1, 1] <- si^2
+    
+    m0 <- rep(0, p)
+    H <- matrix(c(1, rep(0, (p-1))), 1, p)
+    V <- matrix(s0^2)
+    psi <- .psi(psi)
+    
+    ##  Centering: 
+    x <- x - gm$mu
+    ACMres <- ACMfilter(Y=matrix(x,1,N), a=m0, S=Cx, F=Phi, Q=Q, Z=H, V=V, s0=s0, psi=psi, apsi=a, bpsi=b, cpsi=c, flag=flag)
+
+    X.ck <- ACMres$Xf;  X.ck <- X.ck[,2:(N+1)]
+    X   <- ACMres$Xrf; X <- X[,2:(N+1)]
+    st <- as.numeric(unlist(ACMres$rob1L))
+
+
+    if (!lagsmo) {
+        x.ck <- X.ck[1, ]
+        x <- X[1, ]
+    } else {
+        x.ck <- c(X.ck[p, p:N], X.ck[(p-1):1, N])
+        x <- c(X[p, p:N], X[(p-1):1, N])
+    }
+
+##  ARmodel <- .ARmodel(x, p)
+##  y <- ARmodel$y
+##  Z <- ARmodel$Z
+##  r <- resid(lm.fit(Z, y))
+    
+    return(list(filt.ck=x.ck +gm$mu, filt=x + gm$mu, st=st)) #, 
+##              r=c(rep(NA, p), r)))
+
+}

Added: branches/robKalman_bs/pkg/robKalman/R/ACMfilter.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/ACMfilter.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/ACMfilter.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,92 @@
+.getcorrCovACM  <- function (S1, K,  Z, W=diag(nrow(Z)))
+{
+###########################################
+##
+##  R-function: .corrCov - computes filtering error covarince matrix
+##              (internal function)
+##  author: Bernhard Spangl
+##  version: 1.0 (2006-05-22)
+##
+###########################################
+
+##  Paramters:
+##  S1 ... prediction error covariance matrix
+##  K ... Kalman gain
+##  W ... weight matrix
+##  Z ... observation matrix
+
+    S1 - K %*% W %*% Z %*% S1
+}
+
+##steps for classical Kalman filter (cK)
+.ACMinitstep <- function(a, S, ...) 
+              {dots <- list(...)
+               if(hasArg("s0")) 
+                    s0<-dots$"s0"
+               else
+                    s0<-NULL       
+               list( x0 = a,  S0 = S, s0 = s0)}
+
+
+.ACMpredstep <- function (x0, S0, F, Q, rob0, s0, ...)  ### S=P F= Phi
+{
+###########################################
+##
+##  R-function: .ACMpredstep - prediction step (internal function)
+##  author: Bernhard Spangl
+##  version: 1.0 (2006-05-22)
+##
+###########################################
+
+##  Paramters:
+##  x0 ... state vector (filter estimate)
+##  F=Phi ... design matrix of state equation
+##  S0 ... filtering error covariance matrix
+##  Q ... covariance matrix of state innovation process
+##  rob0 ... general robust parameter --- here: scale s0 of nominal Gaussain component of additive noise
+    S1 <- .getpredCov(S0, F, Q)
+    return(list(x1 = F %*% x0, S1 = S1, rob1 = sqrt(S1[1, 1] + s0), Ind=1))
+}
+
+.ACMcorrstep <- function (y, x1, S1, Z, V, rob1, dum=NULL, psi, apsi, bpsi, cpsi, flag, ...)
+{
+###########################################
+##
+##  R-function: .ACMcorrstep - correction step (internal function)
+##  author: Bernhard Spangl
+##  version: 1.0 (2006-05-22)
+##
+###########################################
+
+##  Paramters:
+##  y ... univariate time series 
+##  x1 ... state vector (one-step-ahead predictor)
+##  rob1 ... general robust parameter --- here st ... time-dependent scale parameter
+##  S1 ... prediction error covariance matrix 
+##  Z ... observation matrix
+##  dum ... dummy variable for compatibility with ... argument of calling function
+##  V ... covariance matrix of observation noise
+##  psi ... influence function to be used 
+##  a, b, c ... tuning constants for Hampel's psi-function
+##              (defaul: a=b=2.5, c=5.0)
+##  flag ... character, if "weights" (default), use psi(t)/t to calculate 
+##           the weights; if "deriv", use psi'(t)
+   
+    st <- rob1
+
+    K <- .getKG(S1, Z, V)
+
+    rst <- (y - x1[1])/st
+
+    ps <- psi(rst, apsi, bpsi, cpsi)
+    dx <- K * st * ps
+    x0 <- x1 + dx
+
+    ind <- (abs(rst-ps)>10^-8)
+    
+    w <- psi(rst,  apsi, bpsi, cpsi, flag)
+    
+    S0 <- .getcorrCovACM(S1, K,  Z, W = w*diag(rep(1, nrow(Z))))
+
+    return(list(x0 = x0, K = K,  S0 = S0, Ind=ind, rob0=rob1))
+}

Added: branches/robKalman_bs/pkg/robKalman/R/AllGeneric.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/AllGeneric.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/AllGeneric.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,10 @@
+#if(!isGeneric("type")){ 
+#    setGeneric("type", function(object) standardGeneric("type"))
+#}
+#if(!isGeneric("center")){ 
+#    setGeneric("center", function(object) standardGeneric("center"))
+#}
+#if(!isGeneric("center<-")){
+#    setGeneric("center<-", function(object, value) standardGeneric("center<-"))
+#}
+#

Added: branches/robKalman_bs/pkg/robKalman/R/AllGenerics.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/AllGenerics.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/AllGenerics.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,61 @@
+############################################################################
+# Access methods
+############################################################################
+
+if(!isGeneric("name")) 
+    setGeneric("name", function(object) standardGeneric("name"))
+
+if(!isGeneric("getp")) 
+   setGeneric("getp", function(object) standardGeneric("getp"))
+if(!isGeneric("getq")) 
+   setGeneric("getq", function(object) standardGeneric("getq"))
+
+if(!isGeneric("getF")) 
+   setGeneric("getF", function(object,t) standardGeneric("getF"))
+if(!isGeneric("getZ")) 
+   setGeneric("getZ", function(object,t) standardGeneric("getZ"))
+if(!isGeneric("getQ")) 
+   setGeneric("getQ", function(object,t) standardGeneric("getQ"))
+if(!isGeneric("getV")) 
+   setGeneric("getV", function(object,t) standardGeneric("getV"))
+
+
+############################################################################
+# Replacement methods
+############################################################################
+
+if(!isGeneric("name<-")) 
+    setGeneric("name<-", 
+                function(object, value) standardGeneric("name<-"))
+
+############################################################################
+# generics to  "usual"  methods
+############################################################################
+
+
+# general methods
+
+if(!isGeneric("isOldVersion")) 
+   setGeneric("isOldVersion", function(object) standardGeneric("isOldVersion"))
+
+if(!isGeneric("conv2NewVersion")) 
+   setGeneric("conv2NewVersion", 
+               function(object) standardGeneric("conv2NewVersion"))
+### setting gaps
+
+if(!isGeneric("setgaps"))
+   setGeneric("setgaps", function(object, ...) standardGeneric("setgaps"))
+
+#### generics for log, log10, lgamma, gamma
+
+
+if(!isGeneric("log"))
+   setGeneric("log") #, function(x, base) standardGeneric("log"))
+if(!isGeneric("log10"))
+   setGeneric("log10")
+if(!isGeneric("lgamma"))
+   setGeneric("lgamma")
+if(!isGeneric("gamma"))
+   setGeneric("gamma")
+
+

Added: branches/robKalman_bs/pkg/robKalman/R/AllInitialize.R
===================================================================
--- branches/robKalman_bs/pkg/robKalman/R/AllInitialize.R	                        (rev 0)
+++ branches/robKalman_bs/pkg/robKalman/R/AllInitialize.R	2009-06-09 18:21:49 UTC (rev 30)
@@ -0,0 +1,1012 @@
+#### as to whether to use Generating functions or to use initialize methods:
+#### http://tolstoy.newcastle.edu.au/R/e2/devel/07/01/1976.html
+                     
+################################################################################
+## SPACES
+################################################################################
+
+setMethod("initialize", "Reals",
+          function(.Object) {
+            .Object at dimension <-  1
+            .Object at name <- gettext("Real Space")
+            .Object
+          })
+
+
+setMethod("initialize", "Naturals",
+          function(.Object) {
+            .Object at dimension <-  1
+            .Object at name <- gettext("Grid of Naturals")
+            .Object
+          })
+
+
+################################################################################
+## PARAMETERS
+################################################################################
+
+setMethod("initialize", "GeomParameter",
+          function(.Object, prob = .5) {
+            .Deprecated(new = "new(\"NbinomParameter\"(size = 1, prob, name)",
+                        package = "distr", 
+                        msg = gettext(
+"Class 'GeomParameter' is no longer needed and will be replaced by \nclass 'NbinomParameter' soon."                        
+                        ))
+            .Object at prob <- prob
+            .Object at name <- gettext("Parameter of a Geometric distribution")
+            .Object
+          })
+################################################################################
+## DISTRIBUTIONS
+################################################################################
+
+## Class: UnivariateDistribution
+###produces difficulties in coercing...:
+#
+#setMethod("initialize", "UnivariateDistribution",
+#          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+#                    param = NULL, img = new("Reals"),
+#                    .withSim = FALSE, .withArith = FALSE) {
+#            if(is.null(r)) {
+#              stop("You have at least to give the slot r.")
+#              return(invisible())}
+#            ### Attention: no checking!!!
+#            .Object at img <- img
+#            .Object at param <- param
+#            .Object at d <- d
+#            .Object at p <- p
+#            .Object at q <- q
+#           .Object at r <- r
+#            .Object at .withSim <- .withSim
+#            .Object at .withArith <- .withArith
+#            .Object })
+
+## class AbscontDistribution
+setMethod("initialize", "AbscontDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+                   gaps = NULL, param = NULL, img = new("Reals"),
+                   .withSim = FALSE, .withArith = FALSE) {
+
+            ## don't use this if the call is new("AbscontDistribution")
+            LL <- length(sys.calls())
+            if(sys.calls()[[LL-3]] == "new(\"AbscontDistribution\")")
+               {return(.Object)}
+            
+            if(is.null(r))
+               warning("you have to specify slot r at least")
+                          
+            ## TOBEDONE Errorkanal
+            
+            dpq.approx <- 0
+            
+            dfun <- d
+            pfun <- p
+            qfun <- q
+            
+            if(is.null(d)) {
+              .withSim <- TRUE
+              dpq <- RtoDPQ(r)
+              dpq.approx <- 1
+              dfun <- dpq$dfun}
+            
+            if(is.null(p)) {
+              .withSim <- TRUE
+              if(dpq.approx == 0) {dpq <- RtoDPQ(r)}
+              dpq.approx <- 1
+              pfun <- dpq$pfun}
+            
+            if(is.null(q)) {
+              .withSim <- TRUE
+              if(dpq.approx == 0) {dpq <- RtoDPQ(r)}
+              qfun <- dpq$qfun}
+            
+            .Object at img <- img
+            .Object at param <- param
+            .Object at d <- dfun
+            .Object at p <- pfun
+            .Object at q <- qfun
+            .Object at r <- r
+            .Object at gaps <- gaps
+            .Object at .withSim <- .withSim
+            .Object at .withArith <- .withArith
+            .Object })
+
+## class AffLinAbscontDistribution
+setMethod("initialize", "AffLinAbscontDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, gaps = NULL,
+                   a = 1, b = 0, X0 = Norm(), param = NULL, img = new("Reals"),
+                   .withSim = FALSE, .withArith = FALSE) {
+  X <- new("AbscontDistribution", r = r, d = d, p = p, q = q, gaps = gaps, 
+           param = param, img = img, .withSim = .withSim, 
+           .withArith = .withArith)
+  .Object at gaps  <- X at gaps 
+  .Object at img <- X at img
+  .Object at param <- X at param
+  .Object at a <- a
+  .Object at b <- b
+  .Object at X0 <- X0
+  .Object at d <- X at d
+  .Object at p <- X at p
+  .Object at q <- X at q
+  .Object at r <- X at r
+  .Object at .withSim <- .withSim
+  .Object at .withArith <- .withArith
+  .Object
+})
+
+## Class: DiscreteDistribution
+setMethod("initialize", "DiscreteDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+                    support = NULL, param = NULL, img = new("Reals"), 
+                    .withSim = FALSE, .withArith = FALSE) {
+
+            ## don't use this if the call is new("DiscreteDistribution")
+            LL <- length(sys.calls())
+            if(sys.calls()[[LL-3]] == "new(\"DiscreteDistribution\")")
+               {return(.Object)}
+            
+            if(is.null(r))
+               warning("you have to specify slot r at least")
+              
+            if(is.null(support)) 
+               .Object at support <- as.numeric(names(table(r(10^6))))
+            
+            else .Object at support <- support
+
+            len = length(support)
+
+            if(len > 1){
+              if(min(support[2:len] - support[1:(len - 1)]) < 
+                     getdistrOption("DistrResolution"))
+                 stop("grid too narrow --> change DistrResolution")
+            }
+
+            dpq.approx <- 0
+
+            dfun <- d
+            pfun <- p
+            qfun <- q
+
+            if(is.null(d)) {
+              .withSim <- TRUE
+              dpq <- RtoDPQ.d(r)
+              dpq.approx <- 1
+              dfun <- dpq$dfun
+            }
+
+            if(is.null(p)) {
+              .withSim <- TRUE
+              if(dpq.approx==0) dpq <- RtoDPQ.d(r)
+              dpq.approx <- 1
+              pfun <- dpq$pfun
+            }
+
+            if(is.null(q)) {
+              .withSim <- TRUE
+              if(dpq.approx==0) dpq <- RtoDPQ.d(r)
+              qfun <- dpq$qfun
+            }
+
+            .Object at img <- img
+            .Object at param <- param
+            .Object at d <- dfun
+            .Object at p <- pfun
+            .Object at q <- qfun
+            .Object at r <- r
+            .Object at .withSim <- .withSim
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: AffLinDiscreteDistribution
+setMethod("initialize", "AffLinDiscreteDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+                   support = NULL, a = 1, b = 0, X0 = Binom(), param = NULL, 
+                   img = new("Reals"), .withSim = FALSE, .withArith = FALSE) {
+   ## don't use this if the call is new("DiscreteDistribution")
+   LL <- length(sys.calls())
+   if(sys.calls()[[LL-3]] == "new(\"AffLinDiscreteDistribution\")")
+        X <- new("DiscreteDistribution")
+   else X <- new("DiscreteDistribution", r = r, d = d, p = p, q = q, support = support, 
+             param = param, img = img, .withSim = .withSim, 
+            .withArith = .withArith)
+  .Object at support  <- X at support 
+  .Object at img <- X at img
+  .Object at param <- X at param
+  .Object at a <- a
+  .Object at b <- b
+  .Object at X0 <- X0
+  .Object at d <- X at d
+  .Object at p <- X at p
+  .Object at q <- X at q
+  .Object at r <- X at r
+  .Object at .withSim <- .withSim
+  .Object at .withArith <- .withArith
+  .Object
+})
+
+## Class: LatticeDistribution
+setMethod("initialize", "LatticeDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+                    support = NULL, lattice = NULL, param = NULL, 
+                    img = new("Reals"), .withSim = FALSE, .withArith = FALSE) {
+
+
+             LL <- length(sys.calls())
+             if(sys.calls()[[LL-3]] == "new(\"LatticeDistribution\")")
+             D <- new("DiscreteDistribution")
+             else
+             D <- new("DiscreteDistribution", r = r, d = d, p = p, 
+                       q = q, support = support, param = param, img = img, 
+                     .withSim = .withSim, .withArith = .withArith)
+
+            
+             OS  <- D at support 
+
+             if(is.null(lattice))  
+               {  if(! .is.vector.lattice(OS))
+                      stop("Support as given/generated is not a lattice.")
+                  .Object at lattice <- .make.lattice.es.vector(OS)
+             }else{
+                  .Object at lattice <- lattice 
+             }
+
+
+            .Object at support <- OS
+            .Object at img <- D at img
+            .Object at param <- D at param
+            .Object at d <- D at d
+            .Object at p <- D at p
+            .Object at q <- D at q
+            .Object at r <- D at r
+            .Object at .withSim <- .withSim
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: AffLinLatticeDistribution
+setMethod("initialize", "AffLinLatticeDistribution",
+          function(.Object, r = NULL, d = NULL, p = NULL, q = NULL, 
+                   support = NULL, lattice = NULL, a = 1, b = 0, X0 = Binom(), 
+                   param = NULL, img = new("Reals"), .withSim = FALSE, 
+                   .withArith = FALSE) {
+
+   LL <- length(sys.calls())
+   if(sys.calls()[[LL-3]] == "new(\"AffLinLatticeDistribution\")")
+        X <- new("LatticeDistribution")
+   else X <- new("LatticeDistribution", r = r, d = d, p = p, q = q, 
+                  support = support, lattice = lattice, param = param, 
+                  img = img, .withSim = .withSim, 
+                 .withArith = .withArith)
+
+  .Object at support  <- X at support 
+  .Object at lattice <-  X at lattice 
+  .Object at img <- X at img
+  .Object at param <- X at param
+  .Object at a <- a
+  .Object at b <- b
+  .Object at X0 <- X0
+  .Object at d <- X at d
+  .Object at p <- X at p
+  .Object at q <- X at q
+  .Object at r <- X at r
+  .Object at .withSim <- .withSim
+  .Object at .withArith <- .withArith
+  .Object
+})
+
+######### particular discrete distributions
+
+### Class: Dirac distribution
+setMethod("initialize", "Dirac",
+          function(.Object, location = 0, .withArith = FALSE) {
+            .Object at img <- new("Reals")
+            .Object at param <- new("DiracParameter", location = location)
+            .Object at r <- function(n){}
+            .Object at d <- function(x, log = FALSE){} 
+            .Object at p <- function(q, lower.tail = TRUE, log.p = FALSE){} 
+            .Object at q <- function(p, lower.tail = TRUE, log.p = FALSE){}
+            body(.Object at r) <- substitute({ rep(locationSub, n)},
+                                            list(locationSub = location)
+                                          )
+            body(.Object at d) <- substitute(
+                           { y <- rep(locationSub, length(x))
+                             d0 <- mapply(function(x,y) 
+                                          as.numeric(isTRUE(all.equal(x,y))),
+                                          x = x, y = y)
+                             if (log) d0 <- log(d0)
+                             return(d0)
+                           }, list(locationSub = location)
+                                          )
+            body(.Object at p) <- substitute(
+                           {p0 <-as.numeric(q + .Machine$double.eps^.5 >= 
+                                       locationSub)
+                            if (!lower.tail) p0 <- 1-p0
+                            if (log.p) p0 <- log(p0)
+                            return(p0)
+                            },
+                            list(locationSub = location)
+                                          )
+            body(.Object at q) <- substitute( 
+                { if (log.p) p <- exp(p)
+                  if(any((p < 0)|(p > 1))) 
+                     warning("q Method of class Dirac produced NaN's.")
+                  q0 <- ifelse((p < 0)|(p > 1), NaN, locationSub) 
+                  return(q0)
+                },
+                           list(locationSub = location)
+                                          )
+            .Object at support <- location
+            .Object at lattice <- new("Lattice", pivot = location, width = 1, 
+                                    Length = 1)
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: binomial distribution
+setMethod("initialize", "Binom",
+          function(.Object, size = 1, prob = 0.5, .withArith = FALSE) {
+            .Object at img <- new("Naturals")
+            .Object at param <- new("BinomParameter", size = size, prob = prob)
+            .Object at support <- 0:size
+            .Object at r <- function(n){}
+            .Object at d <- function(x, log = FALSE){}
+            .Object at p <- function(q, lower.tail = TRUE, log.p = FALSE){} 
+            .Object at q <- function(p, lower.tail = TRUE, log.p = FALSE){} 
+            body(.Object at r) <- substitute(
+                           { rbinom(n, size = sizeSub, prob = probSub) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
+            body(.Object at d) <- substitute(
+                           { dbinom(x, size = sizeSub, prob = probSub, 
+                                    log = log) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
+            body(.Object at p) <- substitute(
+                           { pbinom(q, size = sizeSub, prob = probSub, 
+                                    lower.tail = lower.tail, log.p = log.p) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
+            body(.Object at q) <- substitute(
+                           { qbinom(p, size = sizeSub, prob = probSub, 
+                                    lower.tail = lower.tail, log.p = log.p) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
+            .Object at support = 0:size
+            .Object at lattice = new("Lattice", pivot = 0, width = 1, 
+                                   Length = size+1)
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: hypergeometric distribution
+setMethod("initialize", "Hyper",
+          function(.Object, m = 1, n = 1, k = 1, .withArith = FALSE) {
+            .Object at img <- new("Naturals")
+            .Object at param <- new("HyperParameter", m = m, n = n, k = k)
+            .Object at support <- 0:k
+            .Object at r <- function(nn){}
+            .Object at d <- function(x, log = FALSE){}
+            .Object at p <- function(q, lower.tail = TRUE, log.p = FALSE){} 
+            .Object at q <- function(p, lower.tail = TRUE, log.p = FALSE){} 
+            body(.Object at r) <- substitute(
+                               { rhyper(nn, m = mSub, n = nSub, k = kSub) },
+                                 list(mSub = m, nSub = n, kSub = k)
+                                         )
+            body(.Object at d) <- substitute(
+                               { dhyper(x, m = mSub, n = nSub, k = kSub, 
+                                        log = log) },
+                                 list(mSub = m, nSub = n, kSub = k)
+                                          )
+            body(.Object at p) <- substitute(
+                               { phyper(q, m = mSub, n = nSub, k = kSub, 
+                                        lower.tail = lower.tail, log.p = log.p) 
+                                        },
+                                 list(mSub = m, nSub = n, kSub = k)
+                                          )
+            body(.Object at q) <- substitute(
+                               { qhyper(p, m = mSub, n = nSub, k = kSub, 
+                                        lower.tail = lower.tail, log.p = log.p) 
+                                        },
+                                 list(mSub = m, nSub = n, kSub = k)
+                                          )
+            .Object at support <-  seq(from = 0, to = min(k,m), by = 1)
+            .Object at lattice <-  new("Lattice", pivot = 0, width = 1,
+                                     Length = min(k,m)+1 )
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: Poisson distribution 
+setMethod("initialize", "Pois",
+          function(.Object, lambda = 1, .withArith=FALSE) {
+            .Object at img <- new("Naturals")
+            .Object at param <- new("PoisParameter", lambda = lambda)
+            .Object at r <- function(n){}
+            .Object at d <- function(x, log = FALSE){}
+            .Object at p <- function(q, lower.tail = TRUE, log.p = FALSE){} 
+            .Object at q <- function(p, lower.tail = TRUE, log.p = FALSE){} 
+            body(.Object at r) <- substitute({ rpois(n, lambda = lambdaSub) },
+                                            list(lambdaSub = lambda)
+                                          )
+            body(.Object at d) <- substitute({ dpois(x, lambda = lambdaSub, 
+                                                  log = log) },
+                                            list(lambdaSub = lambda)
+                                          )
+            body(.Object at p) <- substitute({ ppois(q, lambda = lambdaSub, 
+                                                  lower.tail = lower.tail, 
+                                                  log.p = log.p) },
+                                            list(lambdaSub = lambda)
+                                          )
+            body(.Object at q) <- substitute({ qpois(p, lambda = lambdaSub, 
+                                                  lower.tail = lower.tail, 
+                                                  log.p = log.p) },
+                                            list(lambdaSub = lambda)
+                                          )
+            .Object at support <- seq(from = 0, by = 1, to = 
+                                   qpois(getdistrOption("TruncQuantile"),
+                                         lambda = lambda, lower.tail = FALSE) 
+                                         + 2
+                                   )
+            .Object at lattice <- new("Lattice", pivot = 0, width = 1, 
+                                    Length = Inf)
+            .Object at .withArith <- .withArith
+            .Object
+          })
+
+## Class: negative binomial distribution
+setMethod("initialize", "Nbinom",
+          function(.Object, size = 1, prob = 0.5, .withArith = FALSE) {
+            .Object at img <- new("Naturals")
+            .Object at param <- new("NbinomParameter", size = size, prob = prob)
+            .Object at r <- function(n){}
+            .Object at d <- function(x, log = FALSE){}
+            .Object at p <- function(q, lower.tail = TRUE, log.p = FALSE){} 
+            .Object at q <- function(p, lower.tail = TRUE, log.p = FALSE){} 
+            body(.Object at r) <- substitute(
+                           { rnbinom(n, size = sizeSub, prob = probSub) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
+            body(.Object at d) <- substitute(
+                           { dnbinom(x, size = sizeSub, prob = probSub, 
+                                     log = log) },
+                             list(sizeSub = size, probSub = prob)
+                                         )
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/robkalman -r 30


More information about the Robkalman-commits mailing list