residual.randomise.lme: Randomise and homogenise residuals from an lme model and add...

Description Usage Arguments Author(s) Examples

View source: R/RepSppFunctions_new.R

Description

Functions used to randomise and homogenise residuals of an lme model.

Usage

1
2

Arguments

mods

List of models

resids

Residuals from model

mod

Individual model for homogenising the residuals.

Author(s)

Robert Bagchi Maintainer: Robert Bagchi <[email protected]>

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
##---- Should be DIRECTLY executable !! ----
##-- ==>  Define data, use random,
##--	or do  help(data=index)  for the standard data sets.

## The function is currently defined as
function (mods, resids) 
{
    level1.resid <- lapply(resids, function(x) x[["level1resids"]])
    N <- max(sapply(level1.resid, length))
    indx <- sample(1:N, replace = T)
    level1.resid.r <- lapply(level1.resid, function(x) return(x[indx]))
    Cmat <- sapply(resids, function(x) attr(x, "zmat"), simplify = FALSE)
    level1.resid.raw.r <- mapply(function(mod, res, Clarge) {
        if (is.null(mod)) 
            return(NULL)
        else {
            wts <- getCovariate(mod$modelStruct$varStruct)[order(order(getGroups(mod)))]
            transform1 <- t(chol(diag(sqrt(wts)) %*% Clarge %*% 
                diag(sqrt(wts))))
            level1.res.raw.r <- as.vector(transform1 %*% res)
            return(level1.res.raw.r)
        }
    }, mod = mods, res = level1.resid.r, Clarge = Cmat, SIMPLIFY = F)
    samp <- lapply(mods, function(mod) {
        if (is.null(mod)) 
            return(NULL)
        else {
            re <- ranef(mod)
            return(lapply(re, function(rj) {
                sample(1:NROW(rj), replace = T)
            }))
        }
    })
    samp <- samp[[1]]
    ranef.r <- mapply(function(mod, res, samp) {
        if (is.null(mod)) 
            return(NULL)
        else {
            ranef.res <- res[-which(names(res) == "level1resids")]
            ranef.res.r <- mapply(function(r, ord) {
                rnew <- as.matrix(r[ord, ])
                rownames(rnew) <- rownames(r)
                return(rnew)
            }, r = ranef.res, ord = samp, SIMPLIFY = FALSE)
            ranef.res.new <- mapply(function(j, r, mod) {
                g <- as.character(getGroups(mod, level = j))
                if (is.null(attr(r, "rownames"))) 
                  rownames(r) <- rownames(ranef(mod, level = j))
                return(r[g, ])
            }, j = as.list(1:length(ranef.res)), r = ranef.res.r, 
                MoreArgs = list(mod = mod), SIMPLIFY = FALSE)
        }
    }, mod = mods, res = resids, MoreArgs = list(samp = samp), 
        SIMPLIFY = FALSE)
    resids <- mapply(function(level1, ranef) {
        list(level1.resid.raw.r = level1, ranef.r = ranef)
    }, level1 = level1.resid.raw.r, ranef = ranef.r, SIMPLIFY = F)
    return(resids)
  }

robertbagchi/ReplicatedPointPatterns documentation built on May 25, 2017, 5:19 a.m.