Статистика Крамера

cramer.R

library(data.table)

#' @title Функция для расчёта статистики Крамера
#' @param x Категориальная переменная.
#' @param y Категориальная переменная.
#' @return Статистика Крамера (число в диапазоне от 0 до 1).
cramer <- function(x, y) {
    # На случай сравнения переменной с самой собой
    if (identical(x, y)) return(1.0)
    # Наблюдаемые частоты
    O <- table(x, y)
    # Количество наблюдений
    n <- length(x)
    n2 <- min(dim(O))
    # Ождиаемые частоты
    E <- outer(rowSums(O), colSums(O), "*") / n
    # Статистика хи-квадрат
    chi <- sum((abs(O - E))^2 / E)
    # Статистика Крамера
    sqrt(chi / (n * (n2 - 1)))
}

#' @title Функция для расчёта статистики Крамера
#' @param x Категориальная переменная.
#' @param y Категориальная переменная.
#' @return Статистика Крамера (число в диапазоне от 0 до 1).
cramer2 <- function(x, y) {
    # На случай сравнения переменной с самой собой
    if (identical(x, y)) return(1.0)
    # Проверка аргументов
    DT <- setDT(list(x = x, y = y))
    n <- DT[, .N]
    # Наблюдаемые частоты
    tb <- DT[, .(o = .N), keyby = .(x, y)]
    tb <- tb[CJ(x, y)]
    if (anyNA(tb[, o)) tb[is.na(o), o := 0]
    # Ожидаемые частоты
    tb[, nx := sum(o), by = x]
    tb[, ny := sum(o), by = y]
    tb[, e := as.numeric(nx) * as.numeric(ny) / sum(n)]
    # Статистика хи-квадрат
    chi <- tb[, sum((abs(o - e))^2 / e)]
    n2 <- tb[, min(sapply(.(x, y), uniqueN))]
    # Статистика Крамера
    sqrt(chi / (n * (n2 - 1)))
}