2015-01-26 2 views
3

R позволяет нам вычислять F-тест между двумя населения:F-тест с P-значения в Python

> d1 = c(2.5579227634, 1.7774243136, 2.0025207896, 1.9518876366, 0.0, 4.1984191803, 5.6170403364, 0.0) 
> d2 = c(16.93800333, 23.2837045311, 1.2674791828, 1.0889208427, 1.0447584137, 0.8971380534, 0.0, 0.0) 
> var.test(d1,d2) 

    F test to compare two variances 

data: d1 and d2 
F = 0.0439, num df = 7, denom df = 7, p-value = 0.000523 
alternative hypothesis: true ratio of variances is not equal to 1 
95 percent confidence interval: 
0.008789447 0.219288957 
sample estimates: 
ratio of variances 
     0.04390249 

Примечание там отчеты P-значение также.

Другой пример, R дал это:

> x1 = c(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 68.7169110318) 
> x2 = c(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.1863361211) 
> var.test(x1,x2) 
#p-value = 1.223e-09 

Что эквивалент в Python? Я проверил это documentation, но, похоже, не дает того, что я хочу.

Этот код дает разные Р-значение (особенно пример 2):

import statistics as stats 
import scipy.stats as ss 
def Ftest_pvalue(d1,d2): 
    """docstring for Ftest_pvalue""" 
    df1 = len(d1) - 1 
    df2 = len(d2) - 1 
    F = stats.variance(d1)/stats.variance(d2) 
    single_tailed_pval = ss.f.cdf(F,df1,df2) 
    double_tailed_pval = single_tailed_pval * 2 
    return double_tailed_pval 

Python дал следующим образом:

In [45]: d1 = [2.5579227634, 1.7774243136, 2.0025207896, 1.9518876366, 0.0, 4.1984191803, 5.6170403364, 0.0] 
In [20]: d2 = [16.93800333, 23.2837045311, 1.2674791828, 1.0889208427, 1.0447584137, 0.8971380534, 0.0, 0.0] 
In [64]: x1 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 68.7169110318] 
In [65]: x2 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.1863361211] 

In [69]: Ftest_pvalue(d1,d2) 
Out[69]: 0.00052297887612346176 

In [70]: Ftest_pvalue(x1,x2) 
Out[70]: 1.9999999987772916 
+2

Скупольное значение p, по-видимому, составляет примерно половину R-единицы (что позволяет использовать представление с плавающей точкой и проблемы округления). Это, казалось бы, подразумевает одностороннее испытание с двумя хвостами. – lvc

+1

@ Комментарий lvc точно верен. Если вы посмотрите на документацию для 'var.test', альтернативная гипотеза по умолчанию - это двухсторонний тест, и когда вы вычисляете' cdf', как вы это делаете с Python, это по сути односторонний тест. –

+0

@lvc: похоже, что после того, как я попробовал другой пример (см. Обновление) – pdubois

ответ

2

rpy2 реализация:

import rpy2.robjects as robjects 
def Ftest_pvalue_rpy2(d1,d2): 
    """docstring for Ftest_pvalue_rpy2""" 
    rd1 = (robjects.FloatVector(d1)) 
    rd2 = (robjects.FloatVector(d2)) 
    rvtest = robjects.r['var.test'] 
    return rvtest(rd1,rd2)[2][0] 

С помощью этого результата:

In [4]: x1 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 68.7169110318] 
In [5]: x2 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.1863361211] 
In [6]: Ftest_pvalue_rpy2(x1,x2) 
Out[6]: 1.2227086010341282e-09 
1

Я должен упомянуть, что xalglib является полный статистических методов, позволяющих сделать это: http://www.alglib.net/ http://www.alglib.net/hypothesistesting/variancetests.php в то время как он является менее гибким, чем оригинальные методы, основанные на SciPy.

Следует отметить, что правильный двойной хвостами процедуру расчета можно найти (в variancetests.c) как:

стата = ae_minreal (xvar/yvar, yvar/xvar, _state); * bothtails = 1- (Распределение Фишера (df1, df2, 1/стат, _state) -fdistribution (df1, df2, стат, _state))

а то, что @Amit Кумар Гупта описывает в своем комментарии ложна (если вы просто удвойте разницу между 1 и односторонним p-значением, вы можете достичь значений выше 1)