Как мне совместно протестировать многоуровневый категориальный эффект в ipython с помощью statsmodels?

Я использую функцию обыкновенных наименьших квадратов (ols) в моделях статистики в ipython, чтобы соответствовать линейной модели, в которой одна ковариата (город) является многоуровневым категориальным эффектом:

result = smf.ols (формула = "Y ~ C (Город) + X * C (Группа)", data = s) .fit ();

(X непрерывно, группа - бинарная категориальная переменная).

Когда я выполняю results.summary (), я получаю по одной строке для каждого уровня города, однако я хотел бы знать общее значение ковариаты «Город» (т. Е. Сравнить Y ~ C (Город) + X * C (Группа) с частичной моделью Y ~ X * C (Группа)).

Есть способ сделать это?

заранее спасибо


person roberto    schedule 05.09.2014    source источник


Ответы (2)


краткий ответ

вы можете использовать anova_lm (тип 3) напрямую или использовать f_test или wald_test и либо построить матрицу ограничений, либо предоставить ограничения гипотезы в виде последовательности формул.

http://statsmodels.sourceforge.net/devel/anova.html http://statsmodels.sourceforge.net/devel/generated/statsmodels.regression.linear_model.RegressionResults.f_test.html

person Josef    schedule 06.09.2014

Спасибо user333700!

Вот подробное описание вашей подсказки. Я генерирую данные с трехуровневой категориальной переменной, использую статистические модели для соответствия модели, а затем совместно тестирую все уровни категориальной переменной:

# 1. generate data
def rnorm(n,u,s):
    return np.random.standard_normal(n)*s+u
a=rnorm(100,-1,1);
b=rnorm(100,0,1);
c=rnorm(100,+1,1);
n=rnorm(300,0,1); # some noise
y=np.concatenate((a,b,c))+n
g=np.zeros(300);
g[0:100]=1
g[100:200]=2
g[200:300]=3
df=pd.DataFrame({'Y':y,'G':g,'N':n});

# 2. fit model
r=smf.ols(formula="Y ~ N + C(G)",data=df).fit();
r.summary()

# 3. joint test
print r.params
A=np.identity(len(r.params)) # identity matrix with size = number of params
GroupTest=A[1:3,:] # for the categorical var., keep the corresponding rows of A
CovTest=A[3,:] # row for the continuous var.
print "Group effect test",r.f_test(GroupTest).fvalue
print "Covariate effect test",r.f_test(CovTest).fvalue

Результат должен быть примерно таким:

Intercept     -1.188975
C(G)[T.2.0]    1.315898
C(G)[T.3.0]    2.137431
N              0.922038
dtype: float64
Group effect test [[ 120.86097747]]
Covariate effect test [[ 259.34155851]]
person roberto    schedule 06.09.2014