Skip to content

单因子CFA ​

检验事先提出的题项与潜变量对应关系。

计算口径 ​

首版为单组、一阶、连续指标的协方差MLW估计。每个因子至少三个题项,各题项仅属于一个因子。采用完整案例;不提供有序题项WLSMV或测量不变性检验。

把同一构念的全部题项选入分析变量。一个因子解释这些指标的协方差,不指定额外误差相关;应至少提供三个指标并核对模型自由度。

结合载荷、CR、AVE与整体拟合判断测量模型。首题载荷固定为1用于定标,不对其显著性作解释。零自由度模型的拟合优度不能证明模型得到支持。

同源实现 ​

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

py
def analyze_structural(data, method, options):
    alpha = number(options, 'alpha', .05, .0001, .25)
    spec = model_description(data, method, options)
    sample = enough(numeric_frame(data, spec['observed']).dropna(), max(20, len(spec['observed'])+2))
    for name in spec['observed']:
        varying(sample[name].values, name)
    correlation = sample.corr().to_numpy()
    if np.min(np.linalg.eigvalsh(correlation)) < 1e-10:
        raise ValueError(_('观测变量相关矩阵奇异,请检查重复指标及完全共线性'))
    renamed = sample.rename(columns=spec['observed_aliases'])
    model = Model(spec['description'])
    free_count = sum(parameter.active for parameter in model.parameters.values())
    moments = len(sample.columns)*(len(sample.columns)+1)//2
    if free_count > moments:
        raise ValueError(_('自由参数超过可用协方差矩,模型欠识别'))
    with warnings.catch_warnings(record=True) as warning_records:
        warnings.simplefilter('always')
        result = model.fit(renamed, obj='MLW', solver='SLSQP', options={'maxiter': 2000, 'ftol': 1e-10})
        if not result.success or not np.all(np.isfinite(model.param_vals)):
            raise ValueError(_('结构方程模型未收敛,请检查量表结构、路径或数据尺度'))
        information = model.calc_fim()
        diagonal = np.diag(information)
        if np.any(diagonal <= 0):
            raise ValueError(_('模型信息矩阵无效,无法计算可靠的标准误'))
        normalized = information/np.sqrt(np.outer(diagonal, diagonal))
        if np.min(np.linalg.eigvalsh(normalized)) <= 1e-9:
            raise ValueError(_('模型局部欠识别或估计不稳定,请减少冗余参数'))
        covariance = np.linalg.inv(information)
        estimates = model.inspect(std_est=True)
        fitted_stats = calc_stats(model).iloc[0].to_dict()
    variances = estimates[(estimates.op == '~~') & (estimates.lval == estimates.rval)]
    if (variances['Estimate'] <= 1e-9).any():
        raise ValueError(_('模型出现接近零或负的误差方差,请检查不当解与题项结构'))
    sigma, ignored_auxiliary = model.calc_sigma()
    # SRMR以样本方差标准化协方差残差,包含唯一协方差元素。
    scale = np.sqrt(np.outer(np.diag(model.mx_cov), np.diag(model.mx_cov)))
    standardized_residuals = (model.mx_cov-sigma)/scale
    lower = np.tril_indices_from(standardized_residuals)
    fitted_stats['SRMR'] = np.sqrt(np.mean(standardized_residuals[lower]**2))
    fit = {key: finite_estimate(value) for key, value in fitted_stats.items()}
    if fit['DoF'] == 0:
        for key in ('RMSEA', 'TLI', 'AGFI', 'chi2 p-value'):
            fit[key] = None
    raw = {'n': len(sample), 'estimator': 'MLW', 'missing_policy': 'listwise', 'model_specification': spec['description'],
           'alias_mapping': spec['inverse'], 'fit': fit, 'estimates': estimates, 'parameter_covariance': covariance,
           'sample_covariance': model.mx_cov, 'implied_covariance': sigma, 'observed_order': model.vars['observed'],
           'warnings': [str(w.message) for w in warning_records], 'free_parameters': free_count}
    measurements, validity, graph_edges = [], [], []
    for factor, items in spec['dimensions']:
        factor_alias = spec['node_aliases'][factor]
        factor_rows = []
        for item in items:
            item_alias = spec['observed_aliases'][item]
            estimate = estimates[(estimates.op == '~') & (estimates.lval == item_alias) & (estimates.rval == factor_alias)].iloc[0]
            loading = float(estimate['Est. Std'])
            if loading*loading >= 1-1e-8:
                raise ValueError(_('标准化因子载荷触及不当解边界,请检查模型'))
            row = {'dimension': factor, 'item': item, 'loading': loading, 'estimate': float(estimate['Estimate']),
                   'se': finite_estimate(estimate['Std. Err']), 'z': finite_estimate(estimate['z-value']),
                   'p': finite_estimate(estimate['p-value'])}
            factor_rows.append(row)
            graph_edges.append({'from': factor_alias, 'to': item_alias, 'label': '%.2f' % loading})
        loadings = np.array([row['loading'] for row in factor_rows])
        errors = 1-loadings**2
        cr = np.sum(loadings)**2/(np.sum(loadings)**2+np.sum(errors))
        ave = np.sum(loadings**2)/(np.sum(loadings**2)+np.sum(errors))
        factor_rows[0].update(cr=cr, ave=ave)
        measurements.extend(factor_rows)
        validity.append({'factor': factor, 'cr': cr, 'ave': ave})
    raw['measurement'] = measurements
    raw['convergent_validity'] = validity
    if spec['dimensions']:
        inner_inverse = np.linalg.inv(np.eye(len(model.mx_beta))-model.mx_beta)
        inner_covariance = inner_inverse @ model.mx_psi @ inner_inverse.T
        aliases = list(model.names_beta[0])
        indexes = [aliases.index(spec['node_aliases'][name]) for name, items in spec['dimensions']]
        latent_covariance = inner_covariance[np.ix_(indexes, indexes)]
        latent_correlation = latent_covariance/np.sqrt(np.outer(np.diag(latent_covariance), np.diag(latent_covariance)))
        discriminant = latent_correlation.copy()
        np.fill_diagonal(discriminant, np.sqrt([row['ave'] for row in validity]))
        raw['factor_correlation'] = latent_correlation
        raw['fornell_larcker'] = {'factors': [name for name, items in spec['dimensions']], 'matrix': discriminant}
    note = _('N=%(n)s;%(fit)s。', n=len(sample), fit=fit_description(fit))
    if method.startswith('stat_cfa_'):
        tables = [PaperTable(_('验证性因子分析'), [Column('dimension', _('维度'), 'text', merge=True), Column('item', _('题项'), 'text'),
                    Column('loading', _('标准化载荷')), Column('estimate', _('非标准化载荷')),
                    Column('se', _('标准误')), Column('z', 'z'), Column('p', 'p', 'p'), Column('cr', 'CR'), Column('ave', 'AVE')], measurements, note)]
        if boolean(options, 'show_discriminant_validity', False):
            factors = [name for name, items in spec['dimensions']]
            rows = [{'factor': name, **{'c'+str(j): discriminant[i, j] for j in range(len(factors))}} for i, name in enumerate(factors)]
            tables.append(PaperTable(_('区分效度'), [Column('factor', _('维度'), 'text')]+[Column('c'+str(j), name) for j, name in enumerate(factors)],
                                     rows, _('对角线为AVE平方根,其余为潜变量相关系数。')))
    else:
        rows = []
        for source, target in spec['edges']:
            estimate = estimates[(estimates.op == '~') & (estimates.lval == spec['node_aliases'][target]) & (estimates.rval == spec['node_aliases'][source])].iloc[0]
            b, se = float(estimate['Estimate']), float(estimate['Std. Err'])
            half = stats.norm.ppf(1-alpha/2)*se
            rows.append({'path': source+' → '+target, 'estimate': b, 'standardized': float(estimate['Est. Std']),
                         'se': se, 'z': float(estimate['z-value']), 'p': float(estimate['p-value']), 'ci': (b-half, b+half)})
            graph_edges.append({'from': spec['node_aliases'][source], 'to': spec['node_aliases'][target], 'label': '%.2f' % float(estimate['Est. Std'])})
        tables = [PaperTable(_('结构路径估计'), [Column('path', _('路径'), 'text'), Column('estimate', 'B'),
                      Column('se', _('标准误')), Column('standardized', 'β'), Column('z', 'z'), Column('p', 'p', 'p'),
                      Column('ci', _('%(level)s%%置信区间', level='%g' % (100*(1-alpha))), 'interval')], rows, note)]
        effects = effect_decomposition(model, spec, estimates, covariance, alpha)
        raw['effects'] = effects
        raw['effect_se_method'] = 'delta_method_full_parameter_covariance'
        if boolean(options, 'show_effects', True) and any(row['effect_code'] == 'indirect' for row in effects):
            tables.append(PaperTable(_('效应分解'), [Column('path', _('路径'), 'text'), Column('effect', _('效应'), 'text'),
                           Column('estimate', _('效应值')), Column('se', _('标准误')), Column('p', 'p', 'p'),
                           Column('ci', _('%(level)s%%置信区间', level='%g' % (100*(1-alpha))), 'interval')], effects,
                           _('效应标准误与区间采用Delta法。')))
        if measurements and boolean(options, 'show_measurement', False):
            tables.append(PaperTable(_('测量模型'), [Column('dimension', _('维度'), 'text'), Column('item', _('题项'), 'text'),
                           Column('loading', _('标准化载荷')), Column('cr', 'CR'), Column('ave', 'AVE')], measurements))
    latent_aliases = set(spec['node_aliases'].values()) if spec['dimensions'] else set()
    image = structural_plot(spec['inverse'], graph_edges, latent_aliases)
    return StatisticalResult(tables, raw, [image])

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

复现本方法 ​

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

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

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

页面操作与解读教程

Released under the AGPL-3.0 License.