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 570, Fri Feb 25 13:21:03 2005 UTC revision 571, Fri Feb 25 13:30:34 2005 UTC
# Line 71  Line 71 
71            }            }
72            fixed.form <- nobars(formula)            fixed.form <- nobars(formula)
73            if (!inherits(fixed.form, "formula")) fixed.form <- ~ 1 # default formula            if (!inherits(fixed.form, "formula")) fixed.form <- ~ 1 # default formula
74              Xmat <- model.matrix(fixed.form, frm)
75            mmats <- c(lapply(random, "[[", 1),            mmats <- c(lapply(random, "[[", 1),
76                       .fixed = list(cbind(model.matrix(fixed.form, frm),                       .fixed = list(cbind(Xmat, .response = model.response(frm))))
77                       .response = model.response(frm))))            ## FIXME: Use Xfrm and Xmat to get the terms and assign
78              ## slots, pass them to lmer_create, then destroy them
79            obj <- .Call("lmer_create", lapply(random, "[[", 2), mmats, PACKAGE = "Matrix")            obj <- .Call("lmer_create", lapply(random, "[[", 2), mmats, PACKAGE = "Matrix")
80              obj@terms <- attr(model.frame(fixed.form, data), "terms")
81              obj@assign <- attr(Xmat, "assign")
82            obj@call <- match.call()            obj@call <- match.call()
83            obj@REML <- REML            obj@REML <- REML
84              rm(Xmat)
85            .Call("lmer_initial", obj, PACKAGE="Matrix")            .Call("lmer_initial", obj, PACKAGE="Matrix")
86            .Call("lmer_ECMEsteps", obj,            .Call("lmer_ECMEsteps", obj,
87                  controlvals$niterEM,                  controlvals$niterEM,
# Line 308  Line 313 
313                            sep = ": "),"")                            sep = ": "),"")
314                return(val)                return(val)
315            } else {            } else {
 #               mCall$terms <- object@terms  
 #               mCall$assign <- object@assign  
316            foo <- object            foo <- object
317            foo@status["factored"] <- FALSE            foo@status["factored"] <- FALSE
318            .Call("lmer_factor", foo, PACKAGE="Matrix")            .Call("lmer_factor", foo, PACKAGE="Matrix")
# Line 319  Line 322 
322            ssr <- ss[[rcol]]            ssr <- ss[[rcol]]
323            ss <- ss[seq(along = dfr)]            ss <- ss[seq(along = dfr)]
324            names(ss) <- object@cnames[[".fixed"]][seq(along = dfr)]            names(ss) <- object@cnames[[".fixed"]][seq(along = dfr)]
325            # FIXME: This only gives single degree of freedom tests                asgn <- foo@assign
326            ms <- ss                terms <- foo@terms
327            df <- rep(1, length(ms))                nmeffects <- attr(terms, "term.labels")
328                  if ("(Intercept)" %in% names(ss))
329                      nmeffects <- c("(Intercept)", nmeffects)
330                  ss <- unlist(lapply(split(ss, asgn), sum))
331                  df <- unlist(lapply(split(asgn,  asgn), length))
332                  dfr <- unlist(lapply(split(dfr, asgn), function(x) x[1]))
333                  ms <- ss/df
334            f <- ms/(ssr/dfr)            f <- ms/(ssr/dfr)
335            P <- pf(f, df, dfr, lower.tail = FALSE)            P <- pf(f, df, dfr, lower.tail = FALSE)
336            table <- data.frame(df, ss, ms, dfr, f, P)            table <- data.frame(df, ss, ms, dfr, f, P)
337            dimnames(table) <-            dimnames(table) <-
338                list(names(ss),                    list(nmeffects,
339                     c("Df", "Sum Sq", "Mean Sq", "Denom", "F value", "Pr(>F)"))                     c("Df", "Sum Sq", "Mean Sq", "Denom", "F value", "Pr(>F)"))
340            if (any(match(names(ss), "(Intercept)", nomatch = 0)))                if ("(Intercept)" %in% nmeffects) table <- table[-1,]
               table <- table[-1,]  
341            attr(table, "heading") <- "Analysis of Variance Table"            attr(table, "heading") <- "Analysis of Variance Table"
342            class(table) <- c("anova", "data.frame")            class(table) <- c("anova", "data.frame")
343            table            table

Legend:
Removed from v.570  
changed lines
  Added in v.571

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