Python中gmodels::fit.contrast的等效实现：药物协同性检验

完全可以用statsmodels实现和R中gmodels::fit.contrast等价的Bliss检验，核心是构造对应对比向量并通过线性模型的对比分析模块计算统计量。以下是具体步骤和代码：

前提假设 #

假设你的TV.data.day数据包含：

• lnTV：肿瘤体积的对数（因变量）
• lnTV.baseline：基线肿瘤体积的对数（协变量）
• Group：分组变量，包含4个水平：对照组（Control）、药物A组（DrugA）、药物B组（DrugB）、联合用药组（Combo）

实现步骤 #

1. 拟合线性模型 #

和R中lm对应，用statsmodels拟合包含协变量和分组的线性模型：

import statsmodels.formula.api as smf

# 拟合模型
model = smf.ols(formula='lnTV ~ lnTV.baseline + Group', data=TV.data.day)
model_fit = model.fit()

2. 构造Bliss检验的对比向量并计算结果 #

Bliss检验的核心是验证联合效应 = 单独效应之和，原假设可转化为：联合组效应 - 药物A效应 - 药物B效应 + 对照组效应 = 0。结合statsmodels默认的处理编码规则，需先确认模型参数顺序，再构造对应对比向量：

# 查看模型参数顺序，确保对比向量对应正确
print(model_fit.params.index)

# 构造对比向量：顺序对应模型参数（截距、基线系数、DrugA哑变量、DrugB哑变量、Combo哑变量）
# 对应Bliss检验的线性组合：Combo效应 - DrugA效应 - DrugB效应
contrast_vec = [0, 0, -1, -1, 1]

# 执行t检验，设置90%置信区间（与R中conf.int=0.9一致）
bliss_test = model_fit.t_test(contrast_vec, alpha=0.1)

# 输出完整检验结果
print("Bliss独立协同性检验结果：")
print(bliss_test.summary())

# 提取关键统计量
print(f"t统计量: {bliss_test.tvalue[0][0]:.4f}")
print(f"p值: {bliss_test.pvalue[0][0]:.4f}")
print(f"90%置信区间: {bliss_test.conf_int(alpha=0.1)[0]}")

结果解读 #

• 若p值小于预设显著性水平（如0.05），则拒绝原假设，说明两种药物存在协同/拮抗作用（联合效应显著不等于单独效应之和）。
• 对比效应值（bliss_test.effect[0][0]）为正代表协同，为负代表拮抗。

关键注意事项 #

• 必须通过model_fit.params.index确认参数顺序，避免因分组编码差异导致对比向量错误。
• statsmodels的t_test默认使用残差自由度，与R中fit.contrast(..., df=TRUE)的统计分布逻辑一致。

内容的提问来源于stack exchange，提问作者Sébastien Wieckowski