Calculate Repeatability Coefficient (from Test-Retest Data)

repeatability.coefficient <- function(x, y)
{
    ## Purpose: Calculate Repeatability Coefficient (from Test-Retest Data)
    ## Arguments:
    ##   x: First Test Outcome (continuous scale)
    ##   y: Second Test Outcome (continuous scale)
    ## Return: Repeatability Coefficient
    ## Author: Feiming Chen
    ## ________________________________________________

    N <- length(x)
    stopifnot(N == length(y))

    s <- sqrt(sum((x - y)^2) / N)
    Repeatability.Coefficient <- 1.96 * s

    a <- paste("Repeatability Coefficient =", round(Repeatability.Coefficient, 3))
    cat("N = ", N, "\n")
    cat("Mean Difference =", round(mean(x-y),3), "\n")
    cat("SD =", round(s, 3), "\n")
    cat(a, "\n")
    plot((x+y)/2, x-y, main = a)
    abline(h=c(0, -Repeatability.Coefficient, Repeatability.Coefficient))
}
if (F) {                                # Unit Test
    N <- 100
    u <- rnorm(N, mean = 100, sd = 5)
    x <-  u + rnorm(N)
    y <-  u + rnorm(N)
    repeatability.coefficient(x, y)     # Theoretical True Value is 2 * sqrt(2) * 1 =  2.8284

    ## Using ANOVA or Regression method for validation
    d <- data.frame(value = c(x, y), group = rep(1:N, 2))
    aov(value ~ factor(group), data = d) # Residual standard error ~ 1
    s <- summary(lm(value ~ factor(group), data = d))$sigma # Residual standard error ~ 1
    ## Y1 - Y2 = E1 - E2, so Var(Y1-Y2) = Var(E1 - E2) = 2 * Var(E1) = 2 * s^2
    cat("Repeatability Coefficient =", 1.96 * sqrt(2) * s, "\n")
}

Example Output:

N =  100 
Mean Difference = 0.012 
SD = 1.433 
Repeatability Coefficient = 2.866