Skip to content

分类变量调节分析 ​

判断X与Y的关系是否随调节变量W而变化。

计算口径 ​

首版包含一个连续X、一个调节变量和连续Y,可加入控制变量。主效应和交互项同时进入模型;连续变量可中心化。

W为分类变量,设置其参照组。交互系数是相对参照组的斜率差;条件效应表分别报告各组X斜率与区间。

先检验交互,再解释条件效应及置信区间。某组显著、另一组不显著不能替代组间斜率差异检验。图中线条为设定条件下的模型关系。

同源实现 ​

以下片段来自 core/statistical/mechanisms.py 的 moderation,由镜像脚本按语法树提取。共享辅助函数和分发逻辑包含在完整下载包中。

py
def moderation(data, method, options):
    predictor = columns(data, options.get('predictors'), 1, 1)[0]
    moderator = columns(data, [options.get('moderator')])[0]
    outcome = options.get('outcome')
    controls = options.get('controls') or []
    if controls:
        columns(data, controls)
    names = [predictor, moderator]+controls
    if len(names) != len(set(names)):
        raise ValueError(_('自变量、调节变量和控制变量不能重复'))
    categorical = options.get('categorical') or []
    if predictor in categorical:
        raise ValueError(_('首版调节分析的焦点自变量须按连续变量处理'))
    is_category = method == 'stat_moderation_categorical'
    if not is_category and moderator in categorical:
        raise ValueError(_('连续调节分析不能将调节变量按分类变量编码'))
    options = {**options, 'categorical': unique_names(categorical+([moderator] if is_category else []))}
    sample = prepare_sample(data, outcome, names, options)
    alpha = number(options, 'alpha', .05, .0001, .25)
    covariance = choice(options, 'covariance', 'standard', ('standard', 'HC1', 'HC3'))
    centered = boolean(options, 'center', True)
    predictor_center = float(sample[predictor].mean()) if centered else 0.
    moderator_center = float(sample[moderator].mean()) if centered and not is_category else 0.
    work = sample.copy()
    work[predictor] -= predictor_center
    if not is_category:
        work[moderator] -= moderator_center
    matrix, terms, references = design_matrix(work, names, options)
    x_index = next(i for i, term in enumerate(terms) if term['variable'] == predictor)
    w_indexes = [i for i, term in enumerate(terms) if term['variable'] == moderator]
    interaction_indexes = []
    for i in w_indexes:
        interaction_indexes.append(matrix.shape[1])
        matrix = np.column_stack([matrix, matrix[:, x_index]*matrix[:, i]])
        terms.append({'key': 'interaction'+str(i), 'label': predictor+' × '+terms[i]['label'],
                      'variable': predictor+' × '+moderator, 'kind': 'interaction'})
    fit = fit_ols(work[outcome].values, matrix, covariance)
    rows = coefficient_rows(fit, matrix, terms, alpha)
    restrictions = np.eye(matrix.shape[1])[interaction_indexes]
    joint = fit.f_test(restrictions)
    if is_category:
        values = [references[moderator]]+[terms[i]['level'] for i in w_indexes]
    else:
        values = options.get('moderator_values')
        if values is None or values == []:
            mean, sd = sample[moderator].mean(), sample[moderator].std(ddof=1)
            values = [mean-sd, mean, mean+sd]
        if not isinstance(values, (list, tuple)) or not 1 <= len(values) <= 8:
            raise ValueError(_('请设置1至8个调节变量取值'))
        try:
            values = [float(v) for v in values]
        except (TypeError, ValueError):
            raise ValueError(_('调节变量取值必须是数字'))
        if not np.all(np.isfinite(values)):
            raise ValueError(_('调节变量取值必须是有限数字'))
    conditional, lines = [], []
    x_grid = np.linspace(sample[predictor].min(), sample[predictor].max(), 45)
    covariance_matrix = fit.cov_params()
    for value in values:
        contrast = np.zeros(matrix.shape[1])
        contrast[x_index] = 1.
        if is_category:
            for index, interaction in zip(w_indexes, interaction_indexes):
                contrast[interaction] = 1. if value == terms[index]['level'] else 0.
        else:
            contrast[interaction_indexes[0]] = value-moderator_center
        estimate = float(contrast @ fit.params)
        variance = float(contrast @ covariance_matrix @ contrast)
        if variance <= 0:
            raise ValueError(_('条件效应的方差无法估计,请检查数据和模型'))
        se = np.sqrt(variance)
        t = estimate/se
        half = stats.t.ppf(1-alpha/2, fit.df_resid)*se
        label = moderator+'='+str(value)
        conditional.append({'term': label, 'b': estimate, 'se': se, 'statistic': t,
                            'p': 2*stats.t.sf(abs(t), fit.df_resid), 'ci': (estimate-half, estimate+half)})
        plot_matrix = np.tile(matrix.mean(axis=0), (len(x_grid), 1))
        plot_matrix[:, 0] = 1
        plot_matrix[:, x_index] = x_grid-predictor_center
        for index, interaction in zip(w_indexes, interaction_indexes):
            plot_matrix[:, index] = (1. if value == terms[index]['level'] else 0.) if is_category else value-moderator_center
            plot_matrix[:, interaction] = plot_matrix[:, x_index]*plot_matrix[:, index]
        lines.append({'label': label, 'x': x_grid, 'y': plot_matrix @ fit.params})
    # 同一张表将回归系数与条件效应并列呈现,复杂诊断留在原始输出。
    for row in conditional:
        rows.append({**row, 'term': _('条件效应:%(value)s', value=row['term'])})
    table = PaperTable(_('调节效应分析'), regression_headers(alpha), rows,
                       _('N=%(n)s;R²=%(r2)s。', n=len(sample), r2='%.3f' % fit.rsquared))
    return StatisticalResult([table], {'model': model_details(fit), 'conditional_effects': conditional,
                                      'interaction_joint_f': float(joint.fvalue), 'interaction_joint_p': float(joint.pvalue),
                                      'interaction_df': len(interaction_indexes), 'references': references,
                                      'centering': {predictor: predictor_center, moderator: moderator_center},
                                      'covariance': covariance, 'parameter_covariance': covariance_matrix},
                             [simple_slopes_plot(lines, predictor, outcome)])

查看完整模块 · 下载算法包 · 复现指南

复现本方法 ​

在解压目录安装 requirements.txt 后执行。样例为固定种子的模拟数据,仅供验证;输出不得冒充真实研究结果。

python
import examples._bootstrap
from examples.statistics_cases import build_data, method_cases
from core.statistical.runner import analyze

method = "stat_moderation_categorical"
options = dict(method_cases())[method]
result = analyze(build_data(), method, options)
print(result.tables[0].html())

页面操作与解读教程

Released under the AGPL-3.0 License.