Разные результаты взвешенного t-теста между R и Python

Я запускаю взвешенный t-тест в Python и вижу разные результаты. Представляется, что проблема заключается в вычислении степени свободы. Хотел бы понять, почему я вижу разные результаты.

Вот пример кода.

In R:

library(weights)
x <- c(373,398,245,272,238,241,134,410,158,125,198,252,577,272,208,260)
y <- c(411,471,320,364,311,390,163,424,228,144,246,371,680,384,279,303)

weightsa = c(rep(1,8), rep(2,8))
weightsb = c(rep(2,8), rep(1,8))

wtd.t.test(x = x,
           y = y,
           weight = weightsa,
           weighty = weightsb, samedata=F)

$test [1] "Two Sample Weighted T-Test (Welch)"

$coefficients
    t.value          df     p.value 
-1.88907197 29.93637837  0.06860382 

$additional Difference     Mean.x     Mean.y   Std. Err   -80.50000 
267.12500  347.62500   42.61352

В Питоне:

import numpy as np
from statsmodels.stats.weightstats import ttest_ind
x = np.asarray([373,398,245,272,238,241,134,410,158,125,198,252,577,272,208,260])
y = np.asarray([411,471,320,364,311,390,163,424,228,144,246,371,680,384,279,303])

weightsa = [1] * 8 + [2] * 8
weightsb = [2] * 8 + [1] * 8

ttest_ind(x, y, usevar='unequal', weights=(weightsa, weightsb))

(-2.3391969704691085, 0.023733058922455107, 45.90244683439944)

Значение P равно 0,06 в R и 0,02 в Python.
В исходном коде R используется формула Саттертуэйта для степеней свободы:

df <- (((vx/n) + (vy/n2))^2)/((((vx/n)^2)/(n - 1)) + 
      ((vy/n2)^2/(n2 - 1)))

Исходный код функции Python также предполагает использование этой формулы:

def dof_satt(self):
        '''degrees of freedom of Satterthwaite for unequal variance
        '''
        d1 = self.d1
        d2 = self.d2
        #this follows blindly the SPSS manual
        #except I use  ``_var`` which has ddof=0
        sem1 = d1._var / (d1.nobs-1)
        sem2 = d2._var / (d2.nobs-1)
        semsum = sem1 + sem2
        z1 = (sem1 / semsum)**2 / (d1.nobs - 1)
        z2 = (sem2 / semsum)**2 / (d2.nobs - 1)
        dof = 1. / (z1 + z2)
        return dof

Числитель здесь выглядит так же, но знаменатель выглядит совершенно иначе.


python r t-test
person matsuo_basho    schedule 24.02.2020    source источник
comment
Знаменатель на самом деле один и тот же, что вы можете подтвердить некоторыми (утомительными) алгебраическими манипуляциями.   -  person duckmayr    schedule 25.02.2020
comment
Так что же тогда отличает вывод?   -  person matsuo_basho    schedule 25.02.2020
comment
Не уверен, но я знаю, что это не степени свободы. Завтра постараюсь продолжить расследование, если к тому времени не будет ответа.   -  person duckmayr    schedule 25.02.2020

Ответы (1)


arrow_upward
3
arrow_downward

Проблема, с которой вы столкнулись здесь, заключается в том, что weights::wtd.t.test() имеет (для меня странный) аргумент по умолчанию mean1 = TRUE, который определяет, «должны ли веса принудительно иметь среднее значение 1» (из help("wtd.t.test")).

Если мы используем mean1 = FALSE, мы получим то же поведение, что и ttest_ind():

wtd.t.test(x = x,
           y = y,
           weight = weightsa,
           weighty = weightsb,
           samedata = FALSE,
           mean1 = FALSE)

$test
[1] "Two Sample Weighted T-Test (Welch)"

$coefficients
    t.value          df     p.value 
-2.33919697 45.90244683  0.02373306 

$additional
Difference     Mean.x     Mean.y   Std. Err 
 -80.50000  267.12500  347.62500   34.41352
person duckmayr    schedule 25.02.2020