SCM

SCM Repository

[matrix] Diff of /pkg/R/lmer.R
ViewVC logotype

Diff of /pkg/R/lmer.R

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 814, Thu Jul 21 23:05:52 2005 UTC revision 820, Wed Jul 27 23:07:04 2005 UTC
# Line 810  Line 810 
810            print(reMat, quote = FALSE)            print(reMat, quote = FALSE)
811        })        })
812    
813  glmmMCMC <- function(obj, nsamp = 1, alpha = 1, beta = 1, saveb = FALSE, verbose = FALSE)  setMethod("mcmcsamp", signature(obj = "lmer"),
814              function(obj, nsamp = 1, verbose = FALSE, saveb = FALSE, ...)
815  {  {
816              if (obj@family$family == "gaussian" &&
817                  obj@family$link == "identity") {
818                  ans <- .Call("lmer_MCMCsamp", obj, saveb, nsamp,
819                                PACKAGE = "Matrix")
820              } else {
821      ## Check arguments      ## Check arguments
     if (!inherits(obj, "lmer")) stop("obj must be of class `lmer'")  
     if (obj@family$family == "gaussian" && obj@family$link == "identity")  
         warn("glmmMCMC not indended for Gaussian family with identity link")  
822      if (length(obj@Omega) > 1 || obj@nc[1] > 1)      if (length(obj@Omega) > 1 || obj@nc[1] > 1)
823          stop("glmmMCMC currently defined for models with a single variance component")                    stop("mcmcsamp currently defined for glmm models with only one variance component")
824      cv <- Matrix:::lmerControl()      cv <- Matrix:::lmerControl()
825      if (verbose) cv$msVerbose <- 1      if (verbose) cv$msVerbose <- 1
826      family <- obj@family      family <- obj@family
# Line 861  Line 864 
864      LMEopt <- getAnywhere("LMEoptimize<-")      LMEopt <- getAnywhere("LMEoptimize<-")
865      doLMEopt <- quote(LMEopt(x = mer, value = cv))      doLMEopt <- quote(LMEopt(x = mer, value = cv))
866      GSpt <- .Call("glmer_init", environment(), PACKAGE = "Matrix")      GSpt <- .Call("glmer_init", environment(), PACKAGE = "Matrix")
     ## FIXME: The definition of shape is not general  
     shape <- length(levels(mer@flist[[1]]))/2 + alpha  
     betainv <- 1/beta  
867      fixed <- obj@fixed      fixed <- obj@fixed
868      varc <- .Call("lmer_coef", mer, 2, PACKAGE = "Matrix")      varc <- .Call("lmer_coef", mer, 2, PACKAGE = "Matrix")
869      b <- .Call("lmer_ranef", mer, PACKAGE = "Matrix")      b <- .Call("lmer_ranef", mer, PACKAGE = "Matrix")
870      ans <- .Call("glmer_MCMCsamp", GSpt, b, fixed, varc, saveb, nsamp,      ans <- .Call("glmer_MCMCsamp", GSpt, b, fixed, varc, saveb, nsamp,
871            PACKAGE = "Matrix")            PACKAGE = "Matrix")
872      .Call("glmer_finalize", GSpt, PACKAGE = "Matrix");      .Call("glmer_finalize", GSpt, PACKAGE = "Matrix");
     class(ans) <- "mcmc"  
     attr(ans, "mcpar") <- c(1, nsamp, 1)  
     ans  
873  }  }
874              cn <- as.character(seq(len = ncol(ans)))
875              fxd <- fixef(obj)
876              cn[seq(along = fxd)] <- names(fxd)
877              cn[length(fxd) + 1] <- "sigma^2"
878              colnames(ans) <- cn
879              ans
880          })
881    
882  rWishart <- function(n, df, invScal)  rWishart <- function(n, df, invScal)
883    .Call("Matrix_rWishart", n, df, invScal)    .Call("Matrix_rWishart", n, df, invScal)
884    

Legend:
Removed from v.814  
changed lines
  Added in v.820

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge