R/Rhat1.R

Defines functions Rhat Rhat1

## File Name: Rhat1.R
## File Version: 1.071


####################################
# auxiliary functions for Rhat statistic
############################################################################
        # Code from rube package
        # Source: http://www.stat.cmu.edu/~hseltman/rube/rube0.2-16/R/Rhat.R

Rhat1 <- function(mat)
{
    m <- ncol(mat)
    n <- nrow(mat)
    b <- apply(mat,2,mean)
    B <- sum((b-mean(mat))^2)*n/(m-1)
    w <- apply(mat,2, stats::var)
    W <- mean(w)
    s2hat <- (n-1)/n*W + B/n
    Vhat <- s2hat + B/m/n
    covWB <- n /m * (stats::cov(w,b^2)-2*mean(b)*stats::cov(w,b))
    varV <- (n-1)^2 / n^2 * stats::var(w)/m +
                (m+1)^2 / m^2 / n^2 * 2*B^2/(m-1) +
                2 * (m-1)*(n-1)/m/n^2 * covWB
    df <- 2 * Vhat^2 / varV
    R <- sqrt((df+3) * Vhat / (df+1) / W)
    return(R)
}



Rhat <- function(arr)
{
    dm <- dim(arr)
    if (length(dm)==2) return(Rhat1(arr))
    if (dm[2]==1) return(NULL)
    if (dm[3]==1) return(Rhat1(arr[,,1]))
    return(apply(arr,3,Rhat1))
}

Try the miceadds package in your browser

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

miceadds documentation built on May 29, 2024, 11:05 a.m.