Использование вычислительных возможностей R для проверки гипотезы о равенстве средних

    Возникла недавно потребность решить вроде бы классическую задачу мат. статистики.
    Проводится испытание определенного push воздействия на группу людей. Необходимо оценить наличие эффекта. Конечно, можно делать это с помощью вероятностного подхода.


    Но рассуждать с бизнесом о нулевых гипотезах и значении p-value совершенно бесполезно и контрпродуктивно.


    Как можно по состоянию на февраль 2019 года сделать это максимально просто и быстро имея под руками ноутбук «средней руки»? Заметка реферативная, формул нет.


    Является продолжением предыдущих публикаций.


    Постановка задачи


    Есть две статистически идентичные по измеряемому показателю группы пользователей (A и B). На группу B оказывается воздействие. Приводит ли это воздействие к изменению среднего значения измеряемому показателю?


    Наиболее популярный вариант — посчитать статистические критерии и сделать вывод. Мне нравится пример "Классические методы статистики: критерий хи-квадрат". При этом совершенно неважно как это выполняется, с помощью спец. программ, Excel, R или еще чего-либо.


    Однако, в достоверности получаемых выводов можно очень сильно сомневаться по следующим причинам:


    1. В действительности мат. статистику мало кто понимает от начала и до конца. Всегда надо держать в голове условия при которых можно применять те или иные методы.
    2. Как правило, использование инструментов и трактовка получаемых результатов идет по приципу однократного вычисления и принятия «светофорного» решения. Чем меньше вопросов, тем лучше для всех участников процесса.

    Критика p-value


    Материалов масса, ссылки на наиболее эффектные из найденных:



    Что можно сделать?


    Сейчас у каждого есть компьютер под руками, поэтому метод Монте-Карло спасает ситуацию. От расчетов p-value переходим к расчету доверительных интервалов (confidence interval) для разницы среднего.


    Книг и материалов множество, но в двух словах (resamapling & fitting) очень компактно изложено в докладе Jake Vanderplas — «Statistics for Hackers» — PyCon 2016. Сама презентация.


    Одна из начальных работ по этой теме, включая предложения по графической визуализации, была написана хорошо известным в советское время популяризатором математики Мартином Гарднером: Confidence intervals rather than P values: estimation rather than hypothesis testing. M.J. Gardner and D.G. Altman, Br Med J (Clin Res Ed). 1986 Mar 15; 292(6522): 746–750.


    Как использовать для этой задачи R?


    Чтобы не делать все руками на нижнем уровне, посмотрим на текущее состояние экосистемы. Не так давно на R был переложен весьма удобный пакет dabestr: Data Analysis using Bootstrap-Coupled Estimation.


    Принципы вычислений и анализа результатов, используемых в dabestr в формате шпаргалок описаны здесь:ESTIMATION STATISTICS BETA ANALYZE YOUR DATA WITH EFFECT SIZES.


    Пример R Notebook для «пощупать»:
    ---
    title: "A/B тестирование средствами bootstrap"
    output: 
      html_notebook:
        self_contained: TRUE
    editor_options: 
      chunk_output_type: inline
    ---

    library(tidyverse)
    library(magrittr)
    library(tictoc)
    library(glue)
    library(dabestr)

    Cимуляция


    Создадим логнормальное распределение длительности операций.


    my_rlnorm <- function(n, mean, sd){
      # пересчитываем мат. моменты: https://en.wikipedia.org/wiki/Log-normal_distribution#Arithmetic_moments
      location <- log(mean^2 / sqrt(sd^2 + mean^2))
      shape <- sqrt(log(1 + (sd^2 / mean^2)))
      print(paste("location:", location))
      print(paste("shape:", shape))
      rlnorm(n, location, shape)  
    }
    
    # N пользователей категории (A = Control)
    A_control <- my_rlnorm(n = 10^3, mean = 500, sd = 150) %T>%
      {print(glue("mean = {mean(.)}; sd = {sd(.)}"))}
    
    # N пользователей категории (B = Test)
    B_test <- my_rlnorm(n = 10^3, mean = 525, sd = 150) %T>%
      {print(glue("mean = {mean(.)}; sd = {sd(.)}"))}

    Собираем данные в виде, необходимом для анализа средствами dabestr, и проводим анализ.


    df <- tibble(Control = A_control, Test = B_test) %>%
      gather(key = "group", value = "value")
    
    tic("bootstrapping")
    two_group_unpaired <- df %>%
      dabest(group, value, 
             # The idx below passes "Control" as the control group, 
             # and "Test" as the test group. The mean difference
             # will be computed as mean(Test) - mean(Control).
             idx = c("Control", "Test"), 
             paired = FALSE,
             reps = 5000
             )
    toc()

    Поглядим на результаты


    two_group_unpaired 
    plot(two_group_unpaired)

    ======================================================


    Результат в виде CI


    DABEST (Data Analysis with Bootstrap Estimation) v0.2.0
    =======================================================
    Unpaired mean difference of Test (n=1000) minus Control (n=1000)
     223 [95CI  209; 236]
    
    5000 bootstrap resamples.
    All confidence intervals are bias-corrected and accelerated.

    и картинки
    image
    вполне понятен и удобен для разговора с бизнесом. Всех расчетов было на «выпить чашечку кофе».


    Предыдущая публикация — «Data Science «спецназ» собственными силами».

    Поделиться публикацией

    Комментарии 18

      0
      Использование вычислительных возможностей R для проверки стат. гипотез

      Вычислительные мощности это же то что хатактеризует «железо», при чем тут язык?
        0

        Допустим, что это утверждение верно. Как в рамках этой гипотезы характеризовать, например, следующие кейсы:


        Кейс #1. Мощный компьютер на котором есть только ОС и компилятор ассемблера. Сможет ли пользователь рассчитать элементарную линейную регрессию по данным, предоставленным в excel файле?


        Кейс #2. К чему относятся вычислительные возможности, предоставляемые SQL (structured query language) ?


        Кейс #3. К чему относятся возможности проведения символьных расчетов в Wolfram Mathematica?


        Кейс #4. Времена 286 процессоров. Неважно, был установлен мат. сопроцессор или нет, операции умножения и деления для чисел с плавающей точкой были доступны в C.


        =============
        R — это экосистема "язык + пакеты", исполняемая поверх определенного программно-аппаратного комплекса (железо + ос).

          0
          Фразы «вычислительные возможности» в контексте языка программмирования — я не встречал вообще. На столько что прочитал как «вычислительные мощности»(каюсь).
          Вроде для языка это будет «выразительные возможности». Такая фраза встречается не редко в литературе по ЯП.

          R — это экосистема «язык + пакеты», исполняемая поверх определенного программно-аппаратного комплекса (железо + ос).

          Как и Октава, Матлаб, Питон(скажем в Анаконде(там кстати и R есть)).
          Но в итоге то всё упирается в железо.
        0
        Всех расчетов было на «выпить чашечку кофе».


        >t.test(x=x,y=y) #всех расчетов на зевнуть…

        :o)
          0

          Ну да, в первом абзаце об этом и было сказано.


          Ровно по этой причине и возникает огромное количество "ложных" открытий. Как правильно применять методы, какие условия должны быть соблюдены, что нужно проверить, как превращать результаты в конкретные действия. Это не только про статистику, а вообще про численные методы.


          Кто из присутствующих возьмет на себя смелость сказать, что он досконально знает мат. статистику? Теперь студенты 2-го, 3-го курсов ведущих вузов не могут объяснить, что такое проекция вектора.

          0

          Проверка статистических гипотез — это формулировка гипотезы, определение уровня значимости, вычисление статистики, построение критических областей и вывод. У вас заголовок со статьей не вяжется.


          То, что бустреппинг хорош для повышения надежности статистических тестов — не новость, но тесты все равно остаются тестами.

            0

            Здесь нет новости, здесь один из кейсов.
            В целом, есть типовая практическая задача по проверке гипотезы о возможности случайного отклонения средних показателей группы A и B. Имея под руками компьютер можно для этого использовать Монте-Карло, результаты проще доносить до бизнеса. Все, никакого скрытого подтекста.

              0

              Значит, я неправильно понял, хотя проверка статистических гипотез — вполне себе сформировавшееся словосочетание, которое обозначает конкретную процедуру.

                0

                я могу переименовать, как было бы лучше?

                  0
                  У вас анализируется конкретный случай использования R для визуализации данных и проверки гипотезы о равенстве средних.
          0
          доверительные интервалы для разницы средних (у вас неверно — «для среднего») можно рассчитать и без бутстрапа. Тем более, что диаграммы вообще не сильно их иллюстрируют, эти диаграммы представляют собой разновидность бокс-плотов (ящиков-с-усами), которые есть в MS Excel, начиная с версии 2016. Для бизнеса можно ими ограничиться.

          Еще по теории. Разница средних — это лишь «сырой» показатель эффекта. Лучше рассчитать d Коэна или g Хеджесса, и построить доверительные интервалы для них.

          Вообще рекомендую книги R.Wilcox. У него есть свой пакет в R — WRS2
            0

            Спасибо, исправил.


            1. Эти представления предложил Мартин Гарднер, надеюсь, многие помнят его книги про математические головоломки. Интересно было ознакомиться с первоисточником, ссылка на статью в тексте.
            2. "Ящик с усами" требует понимания понятия перцентиль. Мой опыт показывает, что далеко не все в бизнесе знают, что это такое.
            3. За ссылку на книги спасибо. Жаль нет времени все прочесть, у меня отложено еще несколько десятков книг. Приходится глядеть урывками по актуальной теме дня. Просто применять пакет\функцию без понимания базовых принципов и идей алгоритмов "под капотом" не очень люблю, можно получить неприятные неожиданности.
              0
              Уточню тогда еще один момент. Кто такой Мартин Гарднер, я не знаю, но понятие доверительных интервалов появилось в 1930-х в работах Е.Неймана и Э.Пирсона. В учебники оно начало проникать в 1950-х. Понятно, что время на книги — дефицит, но рискну порекомендовать книгу ученика Неймана — Эрика Леманна (Erich L.Lehmann) — Fisher, Neyman, and the Creation of Classical Statistics. Без зубодробительной математики. 100 с небольшим страниц, читается легко.
              А главный сторонник доверительных интервалов для меня — это G.Cumming. У него есть прикольный бесплатный образовательный проект для Excel — ESCI, крайне наглядный. На ютьюбе есть его ролики.
            0

            Кроме ошибок первого рода (alpha) в бизнесе могут быть важны и ошибки второго рода (beta). И отсюда понятие мощности теста (power = 1 — beta).


            А вообще, я не советую понял, почему отбутстрапленое распределение выбор.статистики более понятно для бизнеса, чем t-статистика, или p-value, которое в простом варианте трактовки является ошибкой первого рода.

              0
              все-таки лучше не смешивать значение p и уровень ошибок первого рода.

              Бизнесу понятнее картинки как представление данных. Но при этом без эконометрики сегодня бизнес не может. Само по себе понятие доверительного интервала мало кто интерпретирует правильно, это да. Впрочем, как и значение p.

              Cumming G. Understanding replication: Confidence intervals, p-values, and what's likely to happen next time
              Fidler F. Should psychology abandon p-values and teach CIs instead? Evidence-based reforms in statistics education
                0

                Правильно. Но манагеру можно впихнуть: "вероятность, что я тут ошибаюсь, делая вывод о тестовой выборке, равна 1/100К." А если начать, "Это вероятность, что статистика критерия больше наблюдаемого значения при истинности нулевой гипотезы" то совещание просто остановиться. Ха ха ))) Опять же в какой сфере. В Align Tech я репортовал p value, power, effect size. Менеджер был phd from MIT. То есть и бизнес то разный.


                А вот картинку можно сгенеририть и без бутстрапа, просто посчитать SE и построить density средневыборочных или их разницы. Будет даж краше выглядеть, каждый манагеру вешает себе медальку, когда узнает нормальную плотность вероятности! Ха ха ха )))

                  0

                  Бизнес разный, это точно.
                  Хорошо тем, у кого менеджеры думают и понимают, а не сидят на KPI.
                  И нарисовать можно все на миллиметровке и МК-52 в руках.

            Только полноправные пользователи могут оставлять комментарии. Войдите, пожалуйста.

            Самое читаемое