如何在python中对实验组使用OLSResults.f_测试

Question

我试图对我的数据中的三个实验组做一个系数相等的F检验

我已经运行了一个回归来评估随机对照试验的结果，该试验包括四组，G1、G2、G3和对照组

现在我需要确定实验组（G1，G2，G3）是相等的

我知道我可以使用Statsmodel的OLSResults.f_测试来实现这一点。但我不清楚如何配置它。该网站提供了一些例子，但我不知道如何翻译：https://www.statsmodels.org/stable/generated/statsmodels.regression.linear_model.OLSResults.f_test.html

这里给出的例子是：

from statsmodels.datasets import longley
from statsmodels.formula.api import ols
dta = longley.load_pandas().data
formula = 'TOTEMP ~ GNPDEFL + GNP + UNEMP + ARMED + POP + YEAR'
results = ols(formula, dta).fit()
hypotheses = '(GNPDEFL = GNP), (UNEMP = 2), (YEAR/1829 = 1)'
f_test = results.f_test(hypotheses)
print(f_test)

我该如何写下以下假设，以便检查我的3个实验组是否不同

hypotheses = '(G1=G2), (G1=G3), (G2=G3)'

Answer 1

我们可以使用iris示例：

from statsmodels.formula.api import ols
import pandas as pd

data = load_iris()
df = pd.read_csv("http://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data",
                 header=None,names=["s_wid","s_len","p_wid","p_len","species"])

df.species.unique()

array(['Iris-setosa', 'Iris-versicolor', 'Iris-virginica'], dtype=object)

物种分为3类，我们可以像您一样拟合模型：

formula = 's_len ~ species'
results = ols(formula, df).fit()

如果我们看看结果：

results.summary()

Dep. Variable:  s_len   R-squared:  0.392
Model:  OLS Adj. R-squared: 0.384
Method: Least Squares   F-statistic:    47.36
Date:   Sat, 23 May 2020    Prob (F-statistic): 1.33e-16
Time:   01:07:39    Log-Likelihood: -49.688
No. Observations:   150 AIC:    105.4
Df Residuals:   147 BIC:    114.4
Df Model:   2       
Covariance Type:    nonrobust

    coef    std err t   P>|t|   [0.025  0.975]
Intercept   3.4180  0.048   70.998  0.000   3.323   3.513
species[T.Iris-versicolor]  -0.6480 0.068   -9.518  0.000   -0.783  -0.513
species[T.Iris-virginica]   -0.4440 0.068   -6.521  0.000   -0.579  -0.309

如果您的模型仅由上述组组成，则F统计（47.36）和p值（1.33e-16）就是您所需要的。此F-test针对仅拦截模型测试此模型

更详细的解释是：该模型以Iris-setosa为参考进行拟合，另外两个物种对萼片长度的影响s_len计算为与Iris-setosa相关的系数。如果我们看一下平均值，就会明白：

df.groupby('species')['s_len'].mean()

Iris-setosa        3.418
Iris-versicolor    2.770
Iris-virginica     2.974

在这种情况下，假设虹膜花色=0和虹膜维吉尼亚=0，因此各组都相等：

hypotheses = '(species[T.Iris-versicolor] = 0), (species[T.Iris-virginica] = 0)'
results.f_test(hypotheses)

<class 'statsmodels.stats.contrast.ContrastResults'>
<F test: F=array([[47.3644614]]), p=1.327916518456957e-16, df_denom=147, df_num=2>

现在，您可以看到这与摘要中提供的F统计完全相同