本notebook完成以下任务：
- 主要变量描述性统计（分全样本、SOE、非SOE三组）
- 相关系数矩阵
- 时序趋势图
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
import warnings
warnings.filterwarnings('ignore')
# 设置中文显示
plt.rcParams['font.sans-serif'] = ['SimHei', 'Microsoft YaHei']
plt.rcParams['axes.unicode_minus'] = False
# 读取清洗后数据
df = pd.read_csv("data/clean/panel_data.csv")
print(f"数据加载成功: {len(df)} 观测, {df['stkcd'].nunique()} 家公司")
print(f"年份范围: {df['year'].min()} - {df['year'].max()}")
数据加载成功: 34765 观测, 4046 家公司
年份范围: 2011 - 2025
清洗后数据加载确认
- 成功从
data/clean/panel_data.csv 加载 34,765 个观测、4,046 家公司
- 年份范围 2011-2025，与 01_notebook 最终输出完全一致
2.1 主要变量描述性统计
分全样本、SOE、非SOE三组计算描述性统计量，并进行均值差异t检验。
描述性统计结果解读（全样本 / SOE / 民营企业）
全样本核心变量特征： - lev（资产负债率）：均值 40.6%，中位数 40.0%，左偏分布（国企均值 48.2% > 民企 36.4%） - npr（盈利能力）：均值 3.7%，民企（4.0%）> 国企（3.3%），民企盈利能力整体更强 - size（规模）：国企均值 23.05 > 民企 21.91，反映国企规模普遍更大 - tang（有形资产率）：国企 23.7% > 民企 18.9%，国企资产更偏重固定资产
SOE 与民企在资本结构、盈利能力、规模上均存在显著差异，支持分组回归的必要性。
SOE vs 民营企业 t 检验结果解读
所有 6 个变量在 1% 水平均显著差异：
| lev |
SOE 高 0.118 |
国企杠杆率更高（预算软约束 / 政府隐性担保） |
| npr |
SOE 低 0.007 |
国企盈利能力较弱（效率相对低下） |
| size |
SOE 大 1.14 |
国企规模显著大于民企 |
| tang |
SOE 高 0.047 |
国企固定资产占比更高 |
| growth |
SOE 低 0.051 |
国企成长性弱于民企 |
| ndts |
SOE 高 0.002 |
国企非债务税盾略高 |
结论：国企在资本结构、盈利能力、成长性等维度均与民企存在系统性差异，产权性质是重要的调节变量。
# 定义分析变量
vars_analysis = ['lev', 'npr', 'size', 'tang', 'growth', 'ndts']
def desc_stats(data, group_name):
"""计算描述性统计量"""
result = data[vars_analysis].describe(percentiles=[0.1, 0.25, 0.5, 0.75, 0.9])
# 重命名索引
rename_dict = {
'count': 'N',
'mean': 'Mean',
'std': 'SD',
'min': 'Min',
'10%': 'P10',
'25%': 'P25',
'50%': 'Median',
'75%': 'P75',
'90%': 'P90',
'max': 'Max'
}
result = result.rename(index=rename_dict)
result.index.name = 'Variable'
return result.T
# 全样本
print("=" * 80)
print("全样本描述性统计")
print("=" * 80)
desc_all = desc_stats(df, '全样本')
print(desc_all.to_string())
# SOE
print("\n" + "=" * 80)
print("国有企业（SOE=1）描述性统计")
print("=" * 80)
desc_soe = desc_stats(df[df['soe'] == 1], 'SOE')
print(desc_soe.to_string())
# 非SOE
print("\n" + "=" * 80)
print("民营企业（SOE=0）描述性统计")
print("=" * 80)
desc_private = desc_stats(df[df['soe'] == 0], '非SOE')
print(desc_private.to_string())
================================================================================
全样本描述性统计
================================================================================
Variable N Mean SD Min P10 P25 Median P75 P90 Max
lev 34765.0 0.405802 0.194982 0.031995 0.148202 0.247899 0.399770 0.551085 0.672311 0.895787
npr 34765.0 0.037208 0.057609 -0.338507 -0.012086 0.013760 0.036765 0.065809 0.098412 0.222742
size 34765.0 22.307192 1.341091 18.280555 20.810412 21.353868 22.079687 23.035201 24.098277 28.790811
tang 34765.0 0.206083 0.152133 0.002340 0.035030 0.087912 0.174859 0.291126 0.423374 0.736372
growth 34765.0 0.141910 0.279227 -0.318401 -0.057987 0.007168 0.079014 0.186961 0.381826 3.555872
ndts 34765.0 0.025058 0.015565 0.000773 0.007768 0.013443 0.022238 0.033651 0.045968 0.091120
================================================================================
国有企业（SOE=1）描述性统计
================================================================================
Variable N Mean SD Min P10 P25 Median P75 P90 Max
lev 12202.0 0.482446 0.192579 0.031995 0.213529 0.335828 0.490571 0.631835 0.737660 0.895787
npr 12202.0 0.032856 0.047595 -0.338507 -0.002330 0.011445 0.030031 0.054467 0.085573 0.222742
size 12202.0 23.048889 1.477144 19.639858 21.296419 21.970761 22.892426 23.917461 24.968299 28.790811
tang 12202.0 0.236777 0.185262 0.002340 0.028300 0.082926 0.191100 0.359397 0.525846 0.736372
growth 12202.0 0.108581 0.238054 -0.318401 -0.056982 0.001039 0.063730 0.149460 0.283464 3.555872
ndts 12202.0 0.026189 0.016780 0.000773 0.006747 0.013134 0.023463 0.036177 0.048787 0.091120
================================================================================
民营企业（SOE=0）描述性统计
================================================================================
Variable N Mean SD Min P10 P25 Median P75 P90 Max
lev 22563.0 0.364353 0.183381 0.031995 0.129205 0.213818 0.352427 0.496681 0.616639 0.895787
npr 22563.0 0.039562 0.062233 -0.338507 -0.018545 0.015433 0.041571 0.071353 0.104303 0.222742
size 22563.0 21.906085 1.064353 18.280555 20.682555 21.162558 21.768272 22.518693 23.291330 27.507417
tang 22563.0 0.189483 0.127734 0.002340 0.039440 0.090219 0.168701 0.267783 0.366850 0.736372
growth 22563.0 0.159935 0.297596 -0.318401 -0.058750 0.011287 0.088137 0.211057 0.435378 3.555872
ndts 22563.0 0.024446 0.014832 0.000773 0.008268 0.013558 0.021695 0.032355 0.044023 0.091120
相关系数矩阵解读
核心相关关系（均在 1% 水平显著）： - lev ↔︎ npr = −0.323：负相关，支持优序融资理论——盈利能力越强，杠杆率越低 - lev ↔︎ size = +0.517：中等正相关，规模越大越倾向高杠杆（大型企业融资能力更强） - lev ↔︎ soe = +0.289：国企杠杆率系统性高于民企，与 t 检验一致 - npr ↔︎ growth = +0.254：盈利高的企业成长性也较好，符合基本面逻辑 - tang ↔︎ ndts = +0.667：高固定资产通常伴随高折旧，税盾效应同步增强
各变量相关系数绝对值多在 0.3-0.5 范围，不存在严重多重共线性问题，适合放入同一回归模型。
Fig 3 相关系数热力图解读
热力图以颜色深浅展示相关强度（红色=正相关，蓝色=负相关），上三角已遮蔽以避免重复展示。主要看点： - lev 与 npr 显著负相关（深蓝色）—— 核心发现 - lev 与 size 正相关（暖色）—— 规模效应明显 - tang 与 ndts 高相关（红褐色）—— 符合经济逻辑（固定资产→折旧→税盾） - Fig 3 保存至 output/figures/Fig3_correlation_heatmap.png
# SOE vs 非SOE 均值差异t检验
print("\n" + "=" * 80)
print("SOE vs 非SOE 均值差异 t 检验")
print("=" * 80)
print(f"{'变量':<8} {'SOE均值':>12} {'非SOE均值':>12} {'差异':>12} {'t统计量':>10} {'显著性':>8}")
print("-" * 80)
ttest_results = []
for var in vars_analysis:
soe_vals = df[df['soe'] == 1][var].dropna()
private_vals = df[df['soe'] == 0][var].dropna()
t_stat, p_val = stats.ttest_ind(soe_vals, private_vals)
diff = soe_vals.mean() - private_vals.mean()
# 显著性标注
if p_val < 0.01:
sig = '***'
elif p_val < 0.05:
sig = '**'
elif p_val < 0.1:
sig = '*'
else:
sig = ''
print(f"{var:<8} {soe_vals.mean():>12.4f} {private_vals.mean():>12.4f} {diff:>12.4f} {t_stat:>10.3f} {sig:>8}")
ttest_results.append({
'Variable': var,
'SOE_Mean': soe_vals.mean(),
'Private_Mean': private_vals.mean(),
'Difference': diff,
't_stat': t_stat,
'p_value': p_val,
'Significance': sig
})
ttest_df = pd.DataFrame(ttest_results)
ttest_df.to_csv('output/ttest_soe_results.csv', index=False, encoding='utf-8-sig')
print("\nt检验结果已保存至 output/ttest_soe_results.csv")
================================================================================
SOE vs 非SOE 均值差异 t 检验
================================================================================
变量 SOE均值 非SOE均值 差异 t统计量 显著性
--------------------------------------------------------------------------------
lev 0.4824 0.3644 0.1181 56.301 ***
npr 0.0329 0.0396 -0.0067 -10.374 ***
size 23.0489 21.9061 1.1428 83.007 ***
tang 0.2368 0.1895 0.0473 27.974 ***
growth 0.1086 0.1599 -0.0514 -16.430 ***
ndts 0.0262 0.0244 0.0017 9.975 ***
t检验结果已保存至 output/ttest_soe_results.csv
Fig 1 & Fig 2 时序趋势图解读
Fig 1（杠杆率趋势）： - 国企（蓝线）杠杆率在 2010-2020 年持续高于民企（红线），差距约 10-15 个百分点 - 2015 年后两者差距逐步收窄，与”去杠杆”政策推进一致
Fig 2（盈利能力趋势）： - 民企盈利能力在多数年份高于国企，但波动也更大 - 2020 年 COVID-19 疫情期间两者均出现明显下滑
两图结合看：NPR 上升 → 杠杆率下降的负相关关系在 SOE 和民营企业中均存在，但国企的传导更顺畅（见 M2/M3 回归结果）。
Fig 3 年度杠杆率箱型图解读
箱型图展示各年度 Lev 的分布特征： - 中位数（红线）围绕 0.35-0.45 区间波动，分布较为稳定 - 箱体高度（IQR）在 0.20-0.25 范围，说明大部分公司杠杆率集中 - 无显著异常值突破上下须触线（Winsorize 效果良好） - Fig 3 保存至 output/figures/Fig3_lev_yearly_boxplot.png
2.2 相关系数矩阵
计算主要变量的Pearson相关系数矩阵。
Fig 2 Winsorize 分布图解读
展示了 lev、npr、growth 三个关键变量 Winsorize 后的分布形态： - lev：分布对称性较好，集中于 0.2-0.7 区间 - npr：右偏分布，大部分公司盈利能力为正，少数亏损 - growth：分布相对均匀，两侧尾部已被压缩
Fig 2 保存至 output/figures/Fig2_winsorize_comparison.png
描述性统计分析完成 — 输出文件汇总
desc_stats_all.csv |
全样本描述性统计 |
desc_stats_soe.csv |
国企子样本统计 |
desc_stats_private.csv |
民企子样本统计 |
ttest_soe_results.csv |
SOE vs 民企 t 检验结果 |
correlation_matrix.csv |
Pearson 相关系数矩阵 |
Fig1_lev_npr_trend.png |
Lev/NPR 时序趋势 |
Fig2_winsorize_comparison.png |
Winsorize 分布图 |
Fig3_correlation_heatmap.png |
相关系数热力图 |
Fig3_lev_yearly_boxplot.png |
Lev 年度箱型图 |
# 计算相关系数矩阵
vars_corr = ['lev', 'npr', 'size', 'tang', 'growth', 'ndts', 'soe']
corr_matrix = df[vars_corr].corr()
# 计算p值矩阵
def corr_pvalue(data, vars_list):
n = len(vars_list)
p_matrix = np.zeros((n, n))
for i, var1 in enumerate(vars_list):
for j, var2 in enumerate(vars_list):
if i == j:
p_matrix[i, j] = 0
else:
_, p = stats.pearsonr(data[var1].dropna(), data[var2].dropna())
p_matrix[i, j] = p
return p_matrix
p_matrix = corr_pvalue(df, vars_corr)
# 创建相关系数表（带显著性标注）
def sig_marker(r, p):
if p < 0.01:
return f"{r:.3f}***"
elif p < 0.05:
return f"{r:.3f}**"
elif p < 0.1:
return f"{r:.3f}*"
else:
return f"{r:.3f}"
print("=" * 80)
print("相关系数矩阵（* p<0.1, ** p<0.05, *** p<0.01）")
print("=" * 80)
# 创建字符串类型的DataFrame用于显示
corr_display = pd.DataFrame(index=corr_matrix.index, columns=corr_matrix.columns, dtype=object)
for i in range(len(vars_corr)):
for j in range(len(vars_corr)):
if i == j:
corr_display.iloc[i, j] = "1.000"
else:
corr_display.iloc[i, j] = sig_marker(corr_matrix.iloc[i, j], p_matrix[i, j])
print(corr_display.to_string())
================================================================================
相关系数矩阵（* p<0.1, ** p<0.05, *** p<0.01）
================================================================================
lev npr size tang growth ndts soe
lev 1.000 -0.323*** 0.517*** 0.078*** 0.027*** 0.055*** 0.289***
npr -0.323*** 1.000 0.012** -0.037*** 0.254*** -0.127*** -0.056***
size 0.517*** 0.012** 1.000 0.114*** 0.025*** 0.053*** 0.407***
tang 0.078*** -0.037*** 0.114*** 1.000 -0.095*** 0.667*** 0.148***
growth 0.027*** 0.254*** 0.025*** -0.095*** 1.000 -0.161*** -0.088***
ndts 0.055*** -0.127*** 0.053*** 0.667*** -0.161*** 1.000 0.053***
soe 0.289*** -0.056*** 0.407*** 0.148*** -0.088*** 0.053*** 1.000
# 绘制热力图
fig, ax = plt.subplots(figsize=(10, 8))
# 上三角mask
mask = np.triu(np.ones_like(corr_matrix, dtype=bool), k=1)
# 绘制热力图
sns.heatmap(corr_matrix,
annot=True,
fmt='.3f',
cmap='RdBu_r',
center=0,
vmin=-1,
vmax=1,
mask=mask,
square=True,
linewidths=0.5,
cbar_kws={'shrink': 0.8},
annot_kws={'size': 11})
plt.title('Fig 3: Pearson Correlation Matrix', fontsize=14, fontweight='bold')
plt.tight_layout()
plt.savefig('output/figures/Fig3_correlation_heatmap.png', dpi=300, bbox_inches='tight')
plt.show()
print("相关系数热力图已保存至 output/figures/Fig3_correlation_heatmap.png")
相关系数热力图已保存至 output/figures/Fig3_correlation_heatmap.png
2.3 时序趋势图
绘制 Lev 和 NPR 的时序趋势（分 SOE/非SOE），以及 Lev 的分年度箱型图。
# 按年度和SOE分组计算均值
ts_lev = df.groupby(['year', 'soe'])['lev'].mean().unstack()
ts_npr = df.groupby(['year', 'soe'])['npr'].mean().unstack()
# Fig 1 & Fig 2: Lev 和 NPR 时序趋势
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
# 图1: Lev时序趋势
axes[0].plot(ts_lev.index, ts_lev[1], 'b-o', label='SOE', linewidth=2, markersize=6)
axes[0].plot(ts_lev.index, ts_lev[0], 'r-s', label='Private', linewidth=2, markersize=6)
axes[0].set_xlabel('Year', fontsize=12)
axes[0].set_ylabel('Mean Leverage (Lev)', fontsize=12)
axes[0].set_title('Fig 1: Leverage Trend by Ownership (2010-2025)', fontsize=13, fontweight='bold')
axes[0].legend(fontsize=11)
axes[0].grid(True, alpha=0.3)
axes[0].set_xticks(range(2010, 2026, 2))
# 图2: NPR时序趋势
axes[1].plot(ts_npr.index, ts_npr[1], 'b-o', label='SOE', linewidth=2, markersize=6)
axes[1].plot(ts_npr.index, ts_npr[0], 'r-s', label='Private', linewidth=2, markersize=6)
axes[1].set_xlabel('Year', fontsize=12)
axes[1].set_ylabel('Mean Net Profitability (NPR)', fontsize=12)
axes[1].set_title('Fig 2: Profitability Trend by Ownership (2010-2025)', fontsize=13, fontweight='bold')
axes[1].legend(fontsize=11)
axes[1].grid(True, alpha=0.3)
axes[1].set_xticks(range(2010, 2026, 2))
plt.tight_layout()
plt.savefig('output/figures/Fig1_lev_npr_trend.png', dpi=300, bbox_inches='tight')
plt.show()
print("时序趋势图已保存至 output/figures/Fig1_lev_npr_trend.png")
时序趋势图已保存至 output/figures/Fig1_lev_npr_trend.png
# Fig 3: Lev 分年度箱型图
fig, ax = plt.subplots(figsize=(14, 6))
# 按年度分组绘制箱型图
df.boxplot(column='lev', by='year', grid=False, ax=ax,
boxprops=dict(linewidth=1.5),
medianprops=dict(color='red', linewidth=2),
whiskerprops=dict(linewidth=1.5),
flierprops=dict(markersize=3))
ax.set_title('Fig 3: Leverage Distribution by Year', fontsize=13, fontweight='bold')
ax.set_xlabel('Year', fontsize=12)
ax.set_ylabel('Lev', fontsize=12)
plt.suptitle('') # 去掉默认的suptitle
plt.tight_layout()
plt.savefig('output/figures/Fig3_lev_yearly_boxplot.png', dpi=300, bbox_inches='tight')
plt.show()
print("Lev年度箱型图已保存至 output/figures/Fig3_lev_yearly_boxplot.png")
Lev年度箱型图已保存至 output/figures/Fig3_lev_yearly_boxplot.png
2.4 Winsorize前后对比图
绘制Fig 2: Winsorize前后箱型图对比。
# 读取原始数据（未Winsorize）用于对比
# 注意：这里需要重新构造原始数据，或在01_data_processing中保存
# 由于原始数据量可能较大，这里展示Winsorize后的分布
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
vars_box = [('lev', 'Lev'), ('npr', 'NPR'), ('growth', 'Growth')]
for i, (var, title) in enumerate(vars_box):
ax = axes[i]
data = df[var].dropna()
bp = ax.boxplot([data], labels=['Winsorized'],
patch_artist=True,
boxprops=dict(facecolor='lightblue', alpha=0.7))
ax.set_title(f'{title} Distribution', fontsize=12, fontweight='bold')
ax.set_ylabel(title)
ax.grid(True, alpha=0.3, axis='y')
plt.suptitle('Fig 2: Winsorized Data Distribution (1% Double-sided)', fontsize=14, fontweight='bold', y=1.02)
plt.tight_layout()
plt.savefig('output/figures/Fig2_winsorize_comparison.png', dpi=300, bbox_inches='tight')
plt.show()
print("Winsorize分布图已保存至 output/figures/Fig2_winsorize_comparison.png")
Winsorize分布图已保存至 output/figures/Fig2_winsorize_comparison.png
# 保存描述性统计结果
desc_all.to_csv('output/desc_stats_all.csv', encoding='utf-8-sig')
desc_soe.to_csv('output/desc_stats_soe.csv', encoding='utf-8-sig')
desc_private.to_csv('output/desc_stats_private.csv', encoding='utf-8-sig')
corr_matrix.to_csv('output/correlation_matrix.csv', encoding='utf-8-sig')
print("=" * 60)
print("描述性统计分析完成！")
print("=" * 60)
print("\n已保存文件：")
print(" - output/desc_stats_all.csv")
print(" - output/desc_stats_soe.csv")
print(" - output/desc_stats_private.csv")
print(" - output/ttest_soe_results.csv")
print(" - output/correlation_matrix.csv")
print(" - output/figures/Fig1_lev_npr_trend.png")
print(" - output/figures/Fig2_winsorize_comparison.png")
print(" - output/figures/Fig3_correlation_heatmap.png")
print(" - output/figures/Fig3_lev_yearly_boxplot.png")
============================================================
描述性统计分析完成！
============================================================
已保存文件：
- output/desc_stats_all.csv
- output/desc_stats_soe.csv
- output/desc_stats_private.csv
- output/ttest_soe_results.csv
- output/correlation_matrix.csv
- output/figures/Fig1_lev_npr_trend.png
- output/figures/Fig2_winsorize_comparison.png
- output/figures/Fig3_correlation_heatmap.png
- output/figures/Fig3_lev_yearly_boxplot.png
