Python 中 GLM 的方差分析测试

Question

我正在尝试获取 GLM 中每个协变量的 F 统计量和 p 值。在 Python 中，我使用 stats mode.formula.api 来执行 GLM。

formula = 'PropNo_Pred ~ Geography + log10BMI + Cat_OpCavity + CatLes_neles + CatRural_urban + \
        CatPred_Control + CatNative_Intro + Midpoint_of_study'

mod1 = smf.glm(formula=formula, data=A2, family=sm.families.Binomial()).fit()
mod1.summary()

之后，我尝试使用 statsmodels.stats 中的方差分析对该模型进行方差分析测试

table1 = anova_lm(mod3)
print table1

但是我收到一条错误消息：“GLMResults”对象没有属性“ssr”

看起来这个 anova_lm 函数只适用于线性模型，python 中有一个模块可以对 GLM 进行 anova 测试吗？

Answer 1

这是我尝试推出自己的。

嵌套模型的 F 统计量定义为：

(D_s - D_b ) / (addtl_parameters * phi_b)

在哪里：

• D_s是小模型的偏差
• D_b是较大（“大）”模型的偏差
• addtl_parameters是模型之间自由度的差异。
• phi_b是较大模型的色散参数的估计'

“统计理论表明，F 统计量遵循 F 分布，分子自由度等于添加参数的数量，分母自由度等于n - p_b，或者记录数减去大模型中参数的数量”。

我们将其翻译成代码：

from scipy import stats

def calculate_nested_f_statistic(small_model, big_model):
    """Given two fitted GLMs, the larger of which contains the parameter space of the smaller, return the F Stat and P value corresponding to the larger model adding explanatory power"""
    addtl_params = big_model.df_model - small_model.df_model
    f_stat = (small_model.deviance - big_model.deviance) / (addtl_params * big_model.scale)
    df_numerator = addtl_params
    # use fitted values to obtain n_obs from model object:
    df_denom = (big_model.fittedvalues.shape[0] - big_model.df_model)
    p_value = stats.f.sf(f_stat, df_numerator, df_denom)
    return (f_stat, p_value)

这是一个可重现的示例，遵循 statsmodels 中的 gamma GLM 示例（https://www.statsmodels.org/stable/glm.html）：

from scipy import stats

def calculate_nested_f_statistic(small_model, big_model):
    """Given two fitted GLMs, the larger of which contains the parameter space of the smaller, return the F Stat and P value corresponding to the larger model adding explanatory power"""
    addtl_params = big_model.df_model - small_model.df_model
    f_stat = (small_model.deviance - big_model.deviance) / (addtl_params * big_model.scale)
    df_numerator = addtl_params
    # use fitted values to obtain n_obs from model object:
    df_denom = (big_model.fittedvalues.shape[0] - big_model.df_model)
    p_value = stats.f.sf(f_stat, df_numerator, df_denom)
    return (f_stat, p_value)

资料来源： https ://www.casact.org/pubs/monographs/papers/05-Goldburd-Khare-Tevet.pdf

Answer 2

不幸的是，没有。但是，您可以通过对每个项使用模型的假设检验方法来推出自己的模型。事实上，他们的一些ANOVA 方法甚至不使用该属性ssr（这是模型的残差平方和，因此对于二项式 GLM 显然是未定义的）。您可以修改此代码来执行 GLM 方差分析。