1. Уважаемые посетители форума ЭСПП!

    Для просмотра сообщений достаточно прокрутить данное сообщение, а для просмотра списка разделов - вызвать "Каталог".

    Для комментариев необходимо предварительно ознакомиться c Правилами Форума и пройти регистрацию!



    Для того, чтобы быстро ознакомится с возможностями форума, загляните в подраздел Для новичков.

    Если при входе на форум появляется сообщение об ошибке, попробуйте восстановить или сменить пароль, нажав здесь.

Считаем хи-квадрат для таблицы 2х2 в R и на программируемом калькуляторе Casio fx-3650p

Тема в разделе 'Тукачев Ю.А.', создана пользователем Тукачев Ю.А., 28 апр 2016.

  1. Тукачев Ю.А.

    Тукачев Ю.А. Администратор Команда форума

    Исходные данные

    Создадим таблицу 2х2:
    data <-matrix(c(24,13,15,29), nrow = 2)
    Наша таблица 2х2:
    data
    ## [,1] [,2]
    ## [1,] 24 15
    ## [2,] 13 29
    Считаем хи-квадрат с поправкой на непрерывность в R

    chisq.test(data)
    ##
    ## Pearson's Chi-squared test with Yates' continuity correction
    ##
    ## data: data
    ## X-squared = 6.4413, df = 1, p-value = 0.01115
    Считаем хи-квадрат на программируемом калькуляторе

    Программируемый калькулятор Casio fx-3650p имеет память 350 байт для программ. Для вычисления хи-квадрат нам необходимо написать две небольшие программы на языке программирования1:
    • для вычисления эмпирического значения χ2
    • для вычисления p-value (аналог функции pchisq() в R)
    1. Причем язык программирования у калькулятора довольно развитой, в нем присутствуют почти все необходимое для программирования, включая команды для организации циклов и переходов по условию. Программы сохраняются в памяти калькулятора.
    Программа для вычисления эмпирического значения χ2

    Программа занимает 115 байт (115 шагов), а часть ее кода выглядит так:
    115prog.jpg
    Запускаем программу и вводим данные из таблицы 2х2 (A=24, B=15, C=13, D=29). Получаем результат:
    chi.jpg
    Программа для вычисления p-value

    Программа занимает 106 байт (106 шагов), а часть ее кода выглядит так:
    115prog.jpg
    Запускаем программу и вводим данные (A= df = 1, B= эмпир. знач. χ2), получаем результат:
    p_val.jpg
    Видео доступно [по ссылке](https://dl.dropboxusercontent.com/u/33111025/casio/casio.3gp) (2 Мб)
    Мой любимый калькулятор:
    casio.jpg

    Вложения:

    • 106prog.jpg
      106prog.jpg
      Размер файла:
      106,3 КБ
      Просмотров:
      10
    • 115prog.jpg
      115prog.jpg
      Размер файла:
      105,1 КБ
      Просмотров:
      12
    Шмелев А.Г. нравится это.
  2. Шмелев А.Г.

    Шмелев А.Г. Организатор Команда форума

    Юрий, то что R - это актуально, сомневаться не приходится.
    А вот насчет калькулятора я не уверен.


    Вспоминается мне 80-82 гг прошлого века, когда я в практикуме
    на первом курсе учил студентов программировать на калькуляторе
    и даже считать коэффициент ранговой корреляции Спирмена
    на материале Ку-сортировки пятен Роршаха :)


    Сейчас с улыбкой мне вспоминается мой приступ энтузиазма
    тех лет...


    Ваш АШ
  3. Неяскина Ю.Ю.

    Неяскина Ю.Ю. Участник

    Юрий, а, простите, зачем?? Из любви к искусству?))
    Тукачев Ю.А. нравится это.
  4. Тукачев Ю.А.

    Тукачев Ю.А. Администратор Команда форума

    Александр Георгиевич, это я так ностальгирую :) Я учился программировать в конце 80х на микрокалькуляторе Б3-34. Конечно, сейчас мало кто пользуется калькуляторами.
    С другой стороны, удобно же прямо на занятии посчитать хи-квадрат :) или с улыбкой вспомнить "как это было" :)

    Савин Е.Ю. нравится это.
  5. Тукачев Ю.А.

    Тукачев Ю.А. Администратор Команда форума

    Юлия Юрьевна, может быть :)
    Неяскина Ю.Ю. нравится это.
  6. Шмелев А.Г.

    Шмелев А.Г. Организатор Команда форума

    Коллеги, я понял, в чем дело.
    Старая добрая программа на калькуляторе помогает Юрию
    проверить, верно ли работает слишком новая программа на R.
    Так?
  7. Хохлов Н.А.

    Хохлов Н.А. Администратор Команда форума

    Функция chisq.test имеет всё-таки несколько больше возможностей, чем простая формула на программируемом калькуляторе. Например, та же поправка на непрерывность, которая по умолчанию TRUE, но можно ведь и FALSE поставить, и тогда будут совсем иные значения: X-squared = 7.6242, df = 1, p-value = 0.005759. А ещё там есть аргумент simulate.p.value, позволяющий задать симуляцию Монте Карло. В общем, на калькуляторе всего этого не сделать.
    Тукачев Ю.А. нравится это.
  8. Тукачев Ю.А.

    Тукачев Ю.А. Администратор Команда форума

    Я могу поменять форумулу в программе на калькуляторе, чтобы они считал с поправкой или без. Насчет симуляции, то думаю 360 шагов не хватит (после ввода 2х программ свободной памяти 139), чтобы реализовать Монте-Карло (кусок кода в R я выделил жирным), но можно будет попробовать. На самом деле здесь нет разницы между двумя подходами. Если поспмотреть код функции в R:
    function (x, y = NULL, correct = TRUE, p = rep(1/length(x), length(x)),
    rescale.p = FALSE, simulate.p.value = FALSE, B = 2000)
    {
    DNAME <- deparse(substitute(x))
    if (is.data.frame(x))
    x <- as.matrix(x)
    if (is.matrix(x)) {
    if (min(dim(x)) == 1L)
    x <- as.vector(x)
    }
    if (!is.matrix(x) && !is.null(y)) {
    if (length(x) != length(y))
    stop("'x' and 'y' must have the same length")
    DNAME2 <- deparse(substitute(y))
    xname <- if (length(DNAME) > 1L || nchar(DNAME, "w") >
    30)
    ""
    else DNAME
    yname <- if (length(DNAME2) > 1L || nchar(DNAME2, "w") >
    30)
    ""
    else DNAME2
    OK <- complete.cases(x, y)
    x <- factor(x[OK])
    y <- factor(y[OK])
    if ((nlevels(x) < 2L) || (nlevels(y) < 2L))
    stop("'x' and 'y' must have at least 2 levels")
    x <- table(x, y)
    names(dimnames(x)) <- c(xname, yname)
    DNAME <- paste(paste(DNAME, collapse = "\n"), "and",
    paste(DNAME2, collapse = "\n"))
    }
    if (any(x < 0) || any(is.na(x)))
    stop("all entries of 'x' must be nonnegative and finite")
    if ((n <- sum(x)) == 0)
    stop("at least one entry of 'x' must be positive")
    if (simulate.p.value) {
    setMETH <- function() METHOD <<- paste(METHOD, "with simulated p-value\n\t (based on",
    B, "replicates)")
    almost.1 <- 1 - 64 * .Machine$double.eps
    }
    if (is.matrix(x)) {
    METHOD <- "Pearson's Chi-squared test"
    nr <- as.integer(nrow(x))
    nc <- as.integer(ncol(x))
    if (is.na(nr) || is.na(nc) || is.na(nr * nc))
    stop("invalid nrow(x) or ncol(x)", domain = NA)
    sr <- rowSums(x)
    sc <- colSums(x)
    E <- outer(sr, sc, "*")/n
    v <- function(r, c, n) c * r * (n - r) * (n - c)/n^3
    V <- outer(sr, sc, v, n)
    dimnames(E) <- dimnames(x)
    if (simulate.p.value && all(sr > 0) && all(sc > 0)) {
    setMETH()
    tmp <- .Call(C_chisq_sim, sr, sc, B, E)
    STATISTIC <- sum(sort((x - E)^2/E, decreasing = TRUE))
    PARAMETER <- NA
    PVAL <- (1 + sum(tmp >= almost.1 * STATISTIC))/(B +
    1)
    }

    else {
    if (simulate.p.value)
    warning("cannot compute simulated p-value with zero marginals")
    if (correct && nrow(x) == 2L && ncol(x) == 2L) {
    YATES <- min(0.5, abs(x - E))
    if (YATES > 0)
    METHOD <- paste(METHOD, "with Yates' continuity correction")
    }
    else YATES <- 0
    STATISTIC <- sum((abs(x - E) - YATES)^2/E)
    PARAMETER <- (nr - 1L) * (nc - 1L)
    PVAL <- pchisq(STATISTIC, PARAMETER, lower.tail = FALSE)
    }
    }
    else {
    if (length(x) == 1L)
    stop("'x' must at least have 2 elements")
    if (length(x) != length(p))
    stop("'x' and 'p' must have the same number of elements")
    if (any(p < 0))
    stop("probabilities must be non-negative.")
    if (abs(sum(p) - 1) > sqrt(.Machine$double.eps)) {
    if (rescale.p)
    p <- p/sum(p)
    else stop("probabilities must sum to 1.")
    }
    METHOD <- "Chi-squared test for given probabilities"
    E <- n * p
    V <- n * p * (1 - p)
    STATISTIC <- sum((x - E)^2/E)
    names(E) <- names(x)
    if (simulate.p.value) {
    setMETH()
    nx <- length(x)
    sm <- matrix(sample.int(nx, B * n, TRUE, prob = p),
    nrow = n)
    ss <- apply(sm, 2L, function(x, E, k) {
    sum((table(factor(x, levels = 1L:k)) - E)^2/E)
    }, E = E, k = nx)
    PARAMETER <- NA
    PVAL <- (1 + sum(ss >= almost.1 * STATISTIC))/(B +
    1)
    }
    else {
    PARAMETER <- length(x) - 1
    PVAL <- pchisq(STATISTIC, PARAMETER, lower.tail = FALSE)
    }
    }
    names(STATISTIC) <- "X-squared"
    names(PARAMETER) <- "df"
    if (any(E < 5) && is.finite(PARAMETER))
    warning("Chi-squared approximation may be incorrect")
    structure(list(statistic = STATISTIC, parameter = PARAMETER,
    p.value = PVAL, method = METHOD, data.name = DNAME, observed = x,
    expected = E, residuals = (x - E)/sqrt(E), stdres = (x -
    E)/sqrt(V)), class = "htest")
    }

    Последнее редактирование: 28 апр 2016
  9. Тукачев Ю.А.

    Тукачев Ю.А. Администратор Команда форума

    Александр Георгиевич, верно! Плюс осознание того, что я понимаю, что делается с данными, если смог реализовать вычисления на языке калькулятора.