名师讲堂｜使用 R 语言根据上市公司及其子公司与税调匹配结果测算碳排放

2141 字 约 7 分钟 名师讲堂 企业数据R 语言碳排放上市公司数据匹配

指标来源

该碳排放测算方法来自吴超鹏等（2025）发表于《中国工业经济》的论文《高管股权激励与企业碳排放：绿色发展背景下的企业激励机制研究》。该论文采用排放因子法，基于全国税收调查数据（税调）中企业层面的能源消费信息来估算上市公司的碳排放量。

论文原文中对碳排放测算方法的描述（第144页脚注②）如下：

本文的碳排放量指标参考现有文献（Lee et al., 2022; 周泽将等, 2022）的做法，根据全国税收调查数据中的企业用电量、煤炭、石油和天然气的消费量数据，利用排放系数法计算得到。

核心思路是：上市公司自身的能源消费数据通常不可直接获取，但通过将上市公司及其子公司与税调企业进行名称匹配，可以间接获得这些企业的能源消费数据，再利用各能源品种的碳排放系数计算出碳排放量。

计算过程

排放系数

论文中使用的各能源品种碳排放系数如下：

论文中提到的 “每兆瓦时电力的碳排放因子为 7.5 吨” 应该是笔误，10兆瓦时 = 1万kwh，按照《节能低碳技术推广管理暂行办法》里面的描述 0.75kgCO2/kWh = 7.5吨CO2/万千瓦时。

能源品种 碳排放系数 消费量单位
电力 7.5 吨CO2/万千瓦时 万千瓦时
煤炭 1.89 吨CO2/吨 吨
石油 3.02 吨CO2/吨 吨
天然气 21.68 吨CO2/万标立方米 万标立方米

碳排放量的计算公式为：

其中 E 为各能源的消费量，α 为对应的碳排放系数。

数据匹配与汇总流程

整个测算流程分为以下几个步骤：

  1. 读取税调企业能源消费数据（2007~2020年各税调企业电力、煤炭、燃油、天然气消费量面板数据）
  2. 计算各税调企业的碳排放量（消费量 × 碳排放系数）
  3. 对税调企业层面的碳排放进行缩尾处理（1%和99%分位数）
  4. 计算上市公司碳排放（母公司层面）：将税调与上市公司匹配结果合并，按上市公司-年份汇总
  5. 计算全资子公司碳排放：筛选持股比例为100%的全资子公司，与税调匹配后汇总
  6. 合并母公司和子公司碳排放数据，生成最终面板数据

关于全资子公司的筛选

论文中特别指出，对于子公司的碳排放，仅保留持股比例为100%的全资子公司。总持股比例的计算方式为：

总持股比例=直接持股+间接持股

仅保留总持股比例等于100%的子公司，这是为了确保对子公司碳排放的归因准确性。

变量说明

最终生成的数据包含以下变量：

  • total_carbon_emission_ton：总碳排放量（母公司 + 全资子公司），单位：吨CO2
  • ln_carbon_emission：总碳排放量加1取自然对数
  • asinh_carbon_emission：总碳排放量的反双曲正弦变换
  • total_carbon_emission_ton_winsor：缩尾处理后的总碳排放量
  • parent_carbon_emission_ton：母公司层面碳排放量
  • sub_carbon_emission_ton：全资子公司碳排放量
  • 各能源碳排放明细（煤炭、石油、天然气、电力）
  • 各能源消费量（原始单位）

计算代码

下面使用 R 语言完成整个测算过程。需要安装 haven、dplyr、tidyr、labelled 四个包。

加载必要的包

library(haven)
library(dplyr)
library(tidyr)
library(labelled)

第一步：读取数据文件

# 读取税调企业能源消费数据
energy_data <- read_dta("2007～2020年各税调企业电力、煤炭、燃油、天然气消费量面板数据.dta")

# 读取税调与上市公司匹配结果
match_tax_listed <- read_dta("2007～2020年税调与上市公司匹配结果.dta")

# 读取税调与上市公司子公司匹配结果
match_tax_subsidiary <- read_dta("2007~2020年税调数据与上市公司子公司数据匹配结果.dta")

第二步：定义排放系数

# 根据论文中的描述（第144页脚注②）：
# - 每兆瓦时电力的碳排放因子为 7.5 吨 = 7.5 吨CO2/万千瓦时
# - 每吨煤炭的碳排放因子为 1.89 吨 = 1.89 吨CO2/吨煤炭
# - 每吨油的碳排放因子为 3.02 吨 = 3.02 吨CO2/吨油
# - 每万标立方米天然气的碳排放因子为 21.68 吨 = 21.68 吨CO2/万标立方米

ELECTRICITY_EMISSION_FACTOR <- 7.5    # 电力：7.5 吨CO2/万千瓦时
COAL_EMISSION_FACTOR <- 1.89          # 煤炭：1.89 吨CO2/吨
OIL_EMISSION_FACTOR <- 3.02           # 石油：3.02 吨CO2/吨
GAS_EMISSION_FACTOR <- 21.68          # 天然气：21.68 吨CO2/万标立方米

第三步：计算各税调企业的碳排放量

# 此处代码需下载讲义材料查看～

# 查看计算结果
cat("税调企业碳排放统计摘要（吨CO2）：\n")
print(summary(carbon_emission$total_co2_ton))

第四步：缩尾处理

# 此处代码需下载讲义材料查看～

cat("缩尾处理后碳排放统计摘要（吨CO2）：\n")
print(summary(carbon_emission$total_co2_ton_winsor))

第五步：计算上市公司碳排放（母公司层面）

# 将税调企业碳排放数据与上市公司匹配数据合并
# 此处代码需下载讲义材料查看～

cat("母公司层面观测值数量：", nrow(listed_company_carbon), "\n")

第六步：筛选全资子公司并计算碳排放

此处代码需下载讲义材料查看～

第七步：合并母公司和子公司碳排放数据

combined_carbon <- listed_company_carbon %>%
  full_join(
    subsidiary_company_carbon,
    by = c("股票代码" = "证券代码", "年份" = "年份")
  ) %>%
  mutate(
    # 母公司碳排放（如果没有匹配到则为0）
    parent_carbon_emission_ton = ifelse(is.na(carbon_emission_ton), 0, carbon_emission_ton),
    parent_carbon_emission_ton_winsor = ifelse(is.na(carbon_emission_ton_winsor), 0, carbon_emission_ton_winsor),

    # 全资子公司碳排放（如果没有匹配到则为0）
    sub_carbon_emission_ton = ifelse(is.na(subsidiary_carbon_emission_ton), 0, subsidiary_carbon_emission_ton),
    sub_carbon_emission_ton_winsor = ifelse(is.na(subsidiary_carbon_emission_ton_winsor), 0, subsidiary_carbon_emission_ton_winsor),

    # 总碳排放量（母公司 + 全资子公司）
    total_carbon_emission_ton = parent_carbon_emission_ton + sub_carbon_emission_ton,
    total_carbon_emission_ton_winsor = parent_carbon_emission_ton_winsor + sub_carbon_emission_ton_winsor,

    # 取对数后的碳排放量（加1取自然对数）
    ln_carbon_emission = log(total_carbon_emission_ton + 1),
    ln_carbon_emission_winsor = log(total_carbon_emission_ton_winsor + 1),

    # 反双曲正弦变换
    asinh_carbon_emission = asinh(total_carbon_emission_ton),
    asinh_carbon_emission_winsor = asinh(total_carbon_emission_ton_winsor)
  ) %>%
  select(
    股票代码, 股票名称, 年份,
    total_carbon_emission_ton, ln_carbon_emission, asinh_carbon_emission,
    total_carbon_emission_ton_winsor, ln_carbon_emission_winsor, asinh_carbon_emission_winsor,
    parent_carbon_emission_ton, carbon_emission_ton,
    sub_carbon_emission_ton, subsidiary_carbon_emission_ton,
    electricity_consumption, coal_consumption, oil_consumption, gas_consumption,
    sub_elec_cons = subsidiary_electricity_consumption,
    sub_coal_cons = subsidiary_coal_consumption,
    sub_oil_cons = subsidiary_oil_consumption,
    sub_gas_cons = subsidiary_gas_consumption,
    n_matched_firms, n_matched_subsidiaries
  )

第八步：查看结果

cat("样本期间：2007-2020年\n")
cat("上市公司数量：", length(unique(combined_carbon$股票代码)), "\n")
cat("总观测值：", nrow(combined_carbon), "\n")

cat("\n【碳排放量统计摘要】（吨CO2）：\n")
print(summary(combined_carbon$total_carbon_emission_ton))

cat("\n【缩尾处理后碳排放量统计摘要】（吨CO2）：\n")
print(summary(combined_carbon$total_carbon_emission_ton_winsor))

cat("\n【对数碳排放量（ln(Carbon Emission + 1)）统计摘要】：\n")
print(summary(combined_carbon$ln_carbon_emission))

第九步：各年份样本分布

year_dist <- combined_carbon %>%
  group_by(年份) %>%
  summarise(
    公司数量 = n(),
    平均碳排放 = mean(total_carbon_emission_ton_winsor, na.rm = TRUE),
    中位数碳排放 = median(total_carbon_emission_ton_winsor, na.rm = TRUE)
  )
print(year_dist)

第十步：添加变量标签并保存结果

# 为各变量添加中文标签
combined_carbon <- combined_carbon %>%
  labelled::set_variable_labels(
    股票代码        = "上市公司股票代码",
    股票名称        = "上市公司股票名称",
    年份            = "观测年份",
    total_carbon_emission_ton      = "总碳排放量：母公司+全资子公司（吨CO2）",
    ln_carbon_emission             = "总碳排放量加1取自然对数",
    asinh_carbon_emission          = "总碳排放量的反双曲正弦变换",
    total_carbon_emission_ton_winsor = "缩尾处理后的总碳排放量（吨CO2）",
    ln_carbon_emission_winsor      = "缩尾后总碳排放量加1取自然对数",
    asinh_carbon_emission_winsor   = "缩尾后总碳排放量的反双曲正弦变换",
    parent_carbon_emission_ton     = "母公司层面碳排放量（吨CO2）",
    carbon_emission_ton            = "税调匹配的母公司碳排放量（吨CO2）",
    sub_carbon_emission_ton        = "全资子公司碳排放量（吨CO2）",
    subsidiary_carbon_emission_ton = "税调匹配的全资子公司碳排放量（吨CO2）",
    electricity_consumption = "母公司电力消费量（万千瓦时）",
    coal_consumption        = "母公司煤炭消费量（吨）",
    oil_consumption         = "母公司石油消费量（吨）",
    gas_consumption         = "母公司天然气消费量（万标立方米）",
    sub_elec_cons           = "全资子公司电力消费量（万千瓦时）",
    sub_coal_cons           = "全资子公司煤炭消费量（吨）",
    sub_oil_cons            = "全资子公司石油消费量（吨）",
    sub_gas_cons            = "全资子公司天然气消费量（万标立方米）",
    n_matched_firms         = "母公司匹配到的税调企业数量",
    n_matched_subsidiaries  = "匹配到的全资子公司数量"
  )

write_dta(combined_carbon, "2007～2020年上市公司碳排放数据.dta")
cat("数据已保存为 2007～2020年上市公司碳排放数据.dta\n")

点击这里跳转到 RStata 短书平台获取附件：名师讲堂｜使用 R 语言根据上市公司及其子公司与税调匹配结果测算碳排放

相关文章

评论