Files
econometrics/Week7/Week7_Stata代码详解.md

717 lines
44 KiB
Markdown
Raw Permalink Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
# Week 7 · Stata 代码逐段详解
编写与核对:2026-10-26 · 对象:`week7_analysis.do`(165 行) · 工作目录:`C:/Users/cattom/Documents/econometrics/Week7/`
## 目录
- [阅读地图:这份 do 文件完成什么](#overview)
- [语法预备:宏、变量、标量与返回值](#syntax)
- [01 · 初始化环境与记录日志](#setup)
- [02 · 读取 Excel 的指定数值区域](#import)
- [03 · 逐列检查收入表头并核对合计](#audit)
- [04 · 临时切换数据并导出审计表](#export-audit)
- [05 · 从宽表整理成 98 条家庭明细](#reshape)
- [06 · OLS 回归、拟合值与模型统计量](#ols)
- [07 · 手工整理系数、t 检验与置信区间](#coefficient-table)
- [08 · 核对 OLS 残差的两个代数性质](#residual)
- [09 · 按收入水平生成汇总表](#group)
- [10 · 绘制散点与拟合直线,导出四种图形](#plot)
- [11 · 计算异方差稳健标准误并导出](#robust)
- [12 · 恢复 OLS 并进行异方差检验](#hettest)
- [13 · 将模型汇总指标写成一行 CSV](#model-export)
- [14 · 完成标记、关闭日志与退出](#finish)
- [文件索引:每个输出保存了什么](#files)
- [理解与复现时最容易出错的地方](#checkpoints)
- [如何在给定路径下运行和查阅](#run)
- [核对依据与文档范围](#sources)
<a id="overview"></a>
## 阅读地图:这份 do 文件完成什么
这是一份从 Excel 家庭消费明细出发的完整 Stata 分析脚本:核对源表 → 转成长表 → 估计消费与收入的 OLS 关系 → 保存拟合值与残差 → 分组汇总与绘图 → 计算稳健标准误及异方差检验 → 导出结果。回归单位是 98 条家庭记录,收入有 10 个水平。
本文编写和独立核对日期为 **2026 年 10 月 26 日**;所引 Stata 数值来自工作目录内 **2026 年 10 月 23 日**的运行日志及 CSV。此次没有重跑或改动原 do 文件,也没有改写既有分析产物。我们从原 Excel 独立复算了主要统计量,确认它们与现有输出一致。下文代码块直接截取原文件,L 表示原文件物理行号。
> 拟合消费 = 22.930934 + 0.728798742 × 月收入;N = 98;R² = 0.930666
| 输入 / 输出角色 | 内容 |
| --- | --- |
| 原始输入 | `Eg2_1_1/Eg2_1_1.xlsx`,工作表 `Sheet1`,区域 `C3:L17`。 |
| 被解释变量 Y | `consume`:每月消费支出,单位元。 |
| 解释变量 X | `income`:每月可支配收入,单位元。 |
| 主要文件 | DTA 保存数据,CSV 便于通用读取,STER 保存估计结果,GPH/PNG/SVG/PDF 保存图形,LOG 保存运行过程。 |
| 软件与版本 | 用户提供的可执行文件为 `D:/Stata15/StataSE-64.exe`;既有批处理日志显示 Stata/SE 15.1,脚本用 `version 15.0` 指定兼容行为。 |
> **先抓住三个区别**
>
> Excel 行号不等于导入后的 Stata 观测号;`preserve/restore` 管数据,而 `estimates store/restore` 管估计结果;`vce(robust)` 改变标准误和推断,当前模型的 OLS 系数与拟合值保持相同。
<a id="syntax"></a>
## 语法预备:宏、变量、标量与返回值
| 写法 | 在本脚本中的含义 |
| --- | --- |
| `local income = 800` / `` `income' `` | local 定义局部宏;左反引号与右单引号取出宏内容。例如循环首轮 income 宏展开为 800。它不是数据列。 |
| `` `col'[1] `` | 若 col 宏是 C,展开为 C[1]:C 变量的第 1 条观测。方括号是观测下标。 |
| `` consume`income' `` | 拼接变量名:首轮得到 consume800,下一轮得到 consume1100。 |
| `generate income = ...` / `income` | generate 创建数据列;reshape 后的 income 是每条记录都有的数值变量。宏与同名变量用引用语法区分。 |
| `scalar slope_ols = ...` | 保存一个数值;`scalar(slope_ols)` 明确取标量,可避免与变量名混淆。 |
| `r(N)`、`r(sum)`、`r(varlist)` | 一般统计或数据命令返回的结果。summarize 返回 N 和 sum;ds 返回变量名列表。后续 r-class 命令可能覆盖 r(),所以及时保存。 |
| `e(N)`、`e(r2)`、`e(df_r)` | 当前估计命令返回的结果;下一次 regress 会改写当前 e()。 |
| `_b[income]`、`_se[income]`、`_cons` | 当前模型中 income 的系数、标准误,以及截距项的名字。 |
| `_n` 与 `_N` | 分别是当前观测的序号、当前数据集观测总数;删行和排序后含义随当前数据变化。 |
| `=`、`==`、`!`、`^` | 分别常用于赋值/表达式计算、相等比较、逻辑非、乘方。`assert _N == 98` 是检验,不是赋值。 |
| 逗号后的内容 | 一般是命令选项,例如 `regress consume income, vce(robust)`。 |
| `quietly` / `capture` | quietly 隐藏普通输出但仍计算并返回结果;capture 捕获错误,允许程序继续。本文件只在关闭旧日志时使用 capture。 |
| `replace` | 在 save/export/log 等命令中表示覆盖同名文件;它不同于单独的 replace 命令修改变量值。 |
| `tempname` | 为文件句柄、标量或矩阵等对象取得唯一临时名称。本文件把名称放进局部宏,交给 postfile 使用;它不让指定的 DTA 文件变成临时文件。 |
局部宏应连同定义和使用它的代码一起执行;在编辑器中分开执行几段,可能丢失上一段的局部宏上下文。代码中的单引号也不要替换成中文弯引号。语法依据:本机官方帮助 [S1]、[S2]、[S5]。
<a id="setup"></a>
## 01 · 初始化环境与记录日志
**原文件第 1–8 行**
```stata
version 15.0
clear all
set more off
set type double
set linesize 160
capture log close _all
cd "C:/Users/cattom/Documents/econometrics/Week7"
log using "week7_stata.log", text replace
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L1 · `version 15.0` | 按 Stata 15.0 的版本兼容规则解释命令;不是安装或切换可执行软件版本。 |
| L2 · `clear all` | 清理当前内存数据和估计等工作状态,为整段分析建立干净起点。 |
| L3 · `set more off` | 关闭长输出的分页暂停,便于批量运行。 |
| L4 · `set type double` | 未显式指定类型的新数值变量默认使用 double;不等于把所有已有或导入变量自动转成 double。后面的 int 仍按显式类型创建。 |
| L5 · `set linesize 160` | 设置文本输出行宽为 160,减少长表格折行;不改变数据或网页宽度。 |
| L6 · `capture log close _all` | 关闭已打开的所有日志;没有日志时产生的错误由 capture 捕获,不中断程序。 |
| L7 · `cd ...` | 切换工作目录。后续相对路径都以这里为起点;换电脑运行时须修改这个硬编码路径。 |
| L8 · `log using ...` | 把后续命令和结果记录到 week7_stata.log;text 指定纯文本,replace 覆盖同名旧日志。 |
double 是存储类型,format 是显示格式,两者不可互换。增加显示位数不能恢复已经丢失的数值精度。官方帮助:[S1]、[S3]。
<a id="import"></a>
## 02 · 读取 Excel 的指定数值区域
**原文件第 10–12 行**
```stata
display "Week 7: household consumption and disposable income"
display "Input: numeric cells C3:L17; row 17 contains reported totals."
import excel using "Eg2_1_1/Eg2_1_1.xlsx", sheet("Sheet1") cellrange(C3:L17) clear
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L10–11 · `display` | 只打印分析标题和输入范围说明,不生成变量,也不验证数据。 |
| L12 · `import excel using ...` | 从相对路径的工作簿导入 Sheet1 的 C3:L17。clear 允许替换内存数据。 |
| L12 · 没有 `firstrow` | 收入数值表头仍作为第 1 条观测读入;这里变量名为 Excel 列字母 C、D、…、L,后续循环会使用这些名字。 |
| Excel 位置 | 刚导入时的 Stata 观测号 | 含义 |
| --- | --- | --- |
| C3:L3 | 第 1 条 | 收入水平:800、1100、1400、…、3500。 |
| C4:L16 | 第 2–14 条 | 消费明细,13 行 × 10 列;部分为空白。 |
| C17:L17 | 第 15 条 | 原表报告的合计,只用于核对。 |
> **不要把数字表头当成变量名**
>
> 若在这里直接加 firstrow,后续按 C–L 命名和第 1、15 条观测取值的逻辑就需要一起调整。原代码的行号依赖当前导入范围。
导入语法与 firstrow 行为依据本机官方帮助 [S4];具体单元格结构已从当前工作簿核对。
<a id="audit"></a>
## 03 · 逐列检查收入表头并核对合计
**原文件第 14–27 行**
```stata
tempname audit
postfile `audit' double income int n double detail_total reported_total difference using "source_totals_audit.dta", replace
local income = 800
foreach col in C D E F G H I J K L {
assert `col'[1] == `income'
quietly summarize `col' in 2/14
local count = r(N)
local total = r(sum)
local reported = `col'[15]
post `audit' (`income') (`count') (`total') (`reported') (`total' - `reported')
rename `col' consume`income'
local income = `income' + 300
}
postclose `audit'
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L14 · `tempname audit` | 创建唯一的临时句柄名,放在局部宏 audit 中,避免与已有对象重名。 |
| L15 · `postfile ... using ...` | 声明新的审计结果表:income 收入、n 非缺失条数、detail_total 明细和、reported_total 原表合计、difference 差额;n 用 int,其余字段按声明使用 double。不会替换当前 Excel 数据。 |
| L16–17 · `local` + `foreach` | 从收入 800 开始,按 C 到 L 的顺序循环 10 次。每轮 col 宏代表一个列名。 |
| L18 · `assert` | 确认当前列第 1 条观测等于预期收入。若表头或列顺序不符,断言失败,脚本停止。 |
| L19 · `quietly summarize ... in 2/14` | 只汇总明细所在的第 2–14 条观测;缺失值不参加计数和求和;quietly 仅隐藏输出。 |
| L20–22 · 三个 `local` | 保存 r(N)、r(sum) 和第 15 条观测中的原表合计,以免之后命令覆盖 r()。 |
| L23 · `post` | 按声明顺序写入一条审计记录。差额 = 明细之和 − 原表合计;每个待写表达式分别放在括号中。 |
| L24 · `rename` | 例如把 C 改名为 consume800,把 D 改名为 consume1100,为 reshape 的共同前缀 consume 做准备。 |
| L25–26 · 收入递增、结束循环 | 收入宏每轮增加 300,与下一个列的收入表头对应。最后一轮处理 consume3500。 |
| L27 · `postclose` | 完成写入并关闭句柄,使 source_totals_audit.dta 可被后续 use 读入。 |
首轮循环可读成:检查 C[1] = 800 → 统计 C[2] 到 C[14] → 读取 C[15] → 写入审计表 → 把 C 改名为 consume800 → 把收入宏改成 1100。postfile/post/postclose 的行为依据 [S2],summarize 返回值依据 [S5]。
| 收入 | 有效明细数 | 明细之和 | 原表合计 | 差额 |
| --- | --- | --- | --- | --- |
| 800 | 4 | 2420 | 2420 | +0 |
| 1100 | 6 | 4950 | 4950 | +0 |
| 1400 | 11 | 11495 | 11495 | +0 |
| 1700 | 13 | 16445 | 16445 | +0 |
| 2000 | 13 | 19305 | 19305 | +0 |
| 2300 | 13 | 21868 | 23870 | -2002 |
| 2600 | 13 | 25025 | 25025 | +0 |
| 2900 | 10 | 21362 | 21450 | -88 |
| 3200 | 9 | 21153 | 21285 | -132 |
| 3500 | 6 | 15499 | 15510 | -11 |
> **审计记录差异,但不会自动修正**
>
> 有 4 组合计不一致:2300、2900、3200、3500 元组的差额分别为 −2002、−88、−132、−11 元。代码没有 assert difference == 0,因此发现差异后继续分析明细;它既不按合计反推缺失值,也不把合计行加入回归。
<a id="export-audit"></a>
## 04 · 临时切换数据并导出审计表
**原文件第 28–34 行**
```stata
preserve
use "source_totals_audit.dta", clear
list, noobs separator(0)
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "source_totals_audit.csv", datafmt replace
restore
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L28 · `preserve` | 暂存当前已改名的 Excel 宽表。 |
| L29 · `use ... , clear` | 把内存数据切换成刚生成的审计 DTA 文件。 |
| L30 · `list, noobs separator(0)` | 列出审计表;noobs 不显示观测序号,separator(0) 不插入分隔线。 |
| L31 · `quietly ds, has(type numeric)` | 找出数值型变量,将变量名列表放在 r(varlist) 中。 |
| L32 · ``format `r(varlist)' %21.15g`` | 给这些变量设置宽度 21、精度参数 15 的通用数值显示格式 g;必要时使用科学计数法,不是固定 15 位小数。 |
| L33 · `export delimited ... , datafmt replace` | 导出 CSV。datafmt 使用刚设定的显示格式;replace 覆盖同名文件。 |
| L34 · `restore` | 恢复 L28 时保存的宽表;之后继续处理家庭明细。 |
“preserve → 切换/汇总数据 → 导出 → restore”在文件中多次出现,是 Stata 15 单数据集工作流的核心。preserve/restore 不代替 estimates store/restore;数据恢复不应被当作模型恢复。DTA 保存存储值与元数据,CSV 保存文本,显示精度可能影响 CSV 的文本表示。[S3]、[S6]
<a id="reshape"></a>
## 05 · 从宽表整理成 98 条家庭明细
**原文件第 36–53 行**
```stata
keep in 2/14
generate int excel_row = _n + 3
reshape long consume, i(excel_row) j(income)
drop if missing(consume)
sort income excel_row
generate int household_id = _n
order household_id excel_row income consume
assert _N == 98
assert !missing(income, consume)
label variable income "每月可支配收入(元)"
label variable consume "每月消费支出(元)"
label variable excel_row "Source worksheet row"
label variable household_id "Record identifier"
save "households_clean.dta", replace
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "households_clean.csv", datafmt replace
summarize income consume
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L36 · `keep in 2/14` | 仅保留当前宽表第 2–14 条观测:去掉收入表头和合计行,余下 13 行。 |
| L37 · `generate int excel_row = _n + 3` | keep 后 _n 从 1 重新编号;加 3 后得到原 Excel 行号 4–16。这是追溯来源的行标记。 |
| L38 · `reshape long consume, i(excel_row) j(income)` | 将 consume800 至 consume3500 堆叠成长变量 consume;i() 标记原宽表行,j() 把后缀 800…3500 提取成数值变量 income。 |
| L39 · `drop if missing(consume)` | 删除没有消费明细的组合。先形成 13 × 10 = 130 条,再删除 32 条缺失记录,剩 98 条;不会删除数值为 0 的消费。 |
| L40–41 · `sort`、`generate ... household_id` | 按收入、原表行号排序,再以 _n 建立 1–98 的记录编号。编号是人为生成的,不是原始调查家庭 ID。 |
| L42 · `order ...` | 把四个变量按指定顺序排列;不改变数值。 |
| L43 · `assert _N == 98` | 检查总记录数恰好为 98;若数据版本或整理逻辑变化,会停止提醒。 |
| L44 · `assert !missing(income, consume)` | 检查每条观测的收入与消费均非缺失;missing() 在任一参数缺失时为真,! 将其取反。 |
| L45–48 · `label variable` | 设置四个变量的说明标签;变量名仍是 income、consume 等。DTA 保存标签,普通 CSV 不保留完整 Stata 元数据。 |
| L49 · `save ... , replace` | 保存清洗后的 DTA;此时尚未生成拟合值、残差。 |
| L50–52 · `ds`、`format`、`export delimited` | 复用数值变量筛选和高精度 CSV 导出步骤。缩进只是排版,这几行不属于新的代码块。 |
| L53 · `summarize income consume` | 打印收入与消费的个数、均值、标准差、最小值和最大值,进行基本描述统计。 |
| 原 Excel 单元格 | reshape 前 | reshape 后的一条记录 |
| --- | --- | --- |
| C4 = 561 | excel_row=4,consume800=561 | excel_row=4,income=800,consume=561 |
| D4 = 638 | excel_row=4,consume1100=638 | excel_row=4,income=1100,consume=638 |
> **同一 Excel 行不代表同一个家庭**
>
> 第 4 行的不同收入列属于不同家庭消费记录。因此长表中 excel_row 会重复,(excel_row, income) 才标记具体来源单元格。本脚本没有把这批数据设为面板数据,也没有做去重。相同金额可以来自不同家庭。
当前数据平均收入为 2202.0408 元,平均消费为 1627.7755 元。所有 98 条明细逐一与原 Excel 的收入、消费及行号一致。reshape 的 i()/j() 语义依据 [S4]。
<a id="ols"></a>
## 06 · OLS 回归、拟合值与模型统计量
**原文件第 55–69 行**
```stata
regress consume income
estimates store OLS
estimates save "consumption_ols.ster", replace
predict double fitted, xb
predict double residual, residuals
generate double residual_sq = residual^2
scalar slope_ols = _b[income]
scalar intercept_ols = _b[_cons]
scalar rsquared_ols = e(r2)
scalar rmse_ols = e(rmse)
scalar rss_ols = e(rss)
scalar N_ols = e(N)
scalar F_ols = e(F)
scalar df_ols = e(df_r)
scalar pF_ols = Ftail(e(df_m), e(df_r), e(F))
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L55 · `regress consume income` | 第一个变量 consume 是因变量,income 是解释变量;默认含截距,按 98 条明细等权估计 OLS。 |
| L56 · `estimates store OLS` | 把当前估计结果以 OLS 这个名字保存在内存中,供后面恢复。 |
| L57 · `estimates save ...` | 把当前估计结果写入磁盘 STER 文件。它与保存家庭数据的 save 不同,之后可用 estimates use 载入。 |
| L58 · `predict double fitted, xb` | 生成线性预测 fitted = 截距 + 收入系数 × income;xb 表示设计矩阵与系数向量的乘积。 |
| L59 · `predict double residual, residuals` | 生成残差 residual = 实际消费 − fitted。 |
| L60 · `generate double residual_sq = residual^2` | 生成残差平方,便于保留逐户误差信息。后续 estat hettest 自行从估计结果计算所需量,不直接读取这个命名变量。 |
| L61–62 · `_b[income]`、`_b[_cons]` | 分别保存斜率和截距为标量。 |
| L63–66 · `e(r2)`、`e(rmse)`、`e(rss)`、`e(N)` | 分别保存 R²、残差标准差估计(Root MSE)、残差平方和、估计样本量。 |
| L67–68 · `e(F)`、`e(df_r)` | 保存常规模型 F 统计量与残差自由度,本例自由度为 98−2=96。 |
| L69 · `Ftail(...)` | 用 F 分布右尾概率计算整体显著性 p 值;参数依次为模型自由度、残差自由度、F 值,本例为 F(1,96)。 |
> Yᵢ = β₀ + β₁Xᵢ + uᵢ;OLS 最小化 Σ(Yᵢ − β̂₀ − β̂₁Xᵢ)²。
> Ŷᵢ = β̂₀ + β̂₁Xᵢ;êᵢ = Yᵢ − Ŷᵢ;RSS = Σêᵢ²;Root MSE = √[RSS/(N−2)]。
| 统计量 | 现有输出 | 如何理解 |
| --- | --- | --- |
| 截距 / 斜率 | 22.930934 / 0.728798742 | 收入高 100 元,拟合消费高约 72.88 元。0 收入在样本范围外,截距不宜作实质外推。 |
| N / 残差自由度 | 98 / 96 | 估计了截距与一个斜率,共 2 个参数。 |
| R² | 0.930665683 | 收入的线性拟合解释约 93.07% 的消费总离差平方和,不是预测正确率。 |
| Root MSE | 146.627814 元 | 按残差自由度调整后的残差尺度;不是除以 N 的样本内预测 RMSE。 |
| RSS | 2063972.732587 元² | 残差平方和。 |
| F / p_F | 1288.595737 / 1.9519e-57 | 本例整体检验是 H₀: β₁=0,不把截距同时设为 0。 |
拟合值不是实际消费,也不是已知的真实总体条件均值;残差 êᵢ 是估计后可计算的量,真实扰动 uᵢ 通常不可观测。这里的收入系数描述样本中的线性关系;单变量回归本身不建立因果识别。保存标量还可防止后续稳健回归替换 e() 后丢失常规模型统计量。[S7]、[S8]
<a id="coefficient-table"></a>
## 07 · 手工整理系数、t 检验与置信区间
**原文件第 71–84 行**
```stata
tempname coef
postfile `coef' str12 term double estimate double se double t double p double ci_low double ci_high using "regression_coefficients.dta", replace
foreach term in income _cons {
scalar tval = _b[`term'] / _se[`term']
scalar crit = invttail(e(df_r), .025)
post `coef' ("`term'") (_b[`term']) (_se[`term']) (tval) (2*ttail(e(df_r), abs(tval))) (_b[`term']-crit*_se[`term']) (_b[`term']+crit*_se[`term'])
}
postclose `coef'
preserve
use "regression_coefficients.dta", clear
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "regression_coefficients.csv", datafmt replace
restore
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L71–72 · `tempname`、`postfile` | 建立系数结果表。term 用 str12 保存不超过 12 字节的短字符串变量名;其余字段 estimate、se、t、p、ci_low、ci_high 用 double。本例变量名均为 ASCII。 |
| L73 · `foreach term in income _cons` | 循环处理收入系数和截距,最终输出两行结果。 |
| L74 · `tval = _b[...] / _se[...]` | 计算检验系数等于 0 的 t 值;分子没有减去其他假设值,因此原假设是该系数为 0。 |
| L75 · `invttail(e(df_r), .025)` | 得到 t 分布右尾概率 0.025 的临界值,即双侧 95% 区间所用的 t₀.₉₇₅,₉₆;.025 不表示区间置信度为 2.5%。 |
| L76 · `post` | 依次写入变量名、系数、标准误、t 值、双侧 p 值,以及置信区间上下限。abs(tval) 取绝对值,ttail 求右尾概率,再乘 2。 |
| L77–78 · 循环结束、`postclose` | 完成两行写入,关闭系数表。 |
| L79–84 · `preserve` 至 `restore` | 暂存家庭明细,读取系数表、设置数值格式并导出 CSV,再恢复家庭明细。模型与数据是两个层面的对象。 |
> t = β̂ / SE(β̂);p = 2 × P(T₉₆ ≥ |t|);95% CI = β̂ ± t₀.₉₇₅,₉₆ × SE(β̂)。
| 变量 | 估计系数 | 常规标准误 | t | p | 95% 置信区间 |
| --- | --- | --- | --- | --- | --- |
| income | 0.728799 | 0.020302 | 35.8970 | 1.95192e-57 | [0.688499, 0.769099] |
| _cons | 22.930934 | 47.096629 | 0.4869 | 0.627445 | [-70.555136, 116.417004] |
临界值约为 1.984984。收入系数显著为正;截距的常规 p 值约 0.627,不能拒绝截距为 0。很小的 p 值不等于系数具有很大的经济效应。当前采用常规 OLS 标准误;有限样本精确 t/F 推断还需要相应模型条件,不能只看公式就认定条件成立。[S2]、[S7]、[S9]
<a id="residual"></a>
## 08 · 核对 OLS 残差的两个代数性质
**原文件第 86–97 行**
```stata
quietly summarize residual, meanonly
scalar residual_sum = r(sum)
generate double income_residual = income * residual
quietly summarize income_residual, meanonly
scalar income_residual_sum = r(sum)
display "Sum of residuals = " %21.12g scalar(residual_sum)
display "Sum of income * residual = " %21.12g scalar(income_residual_sum)
drop income_residual
save "households_fitted.dta", replace
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "households_fitted.csv", datafmt replace
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L86–87 · `summarize residual, meanonly` | meanonly 省去方差计算并抑制汇总显示,但仍返回 r(sum);把残差和存入标量 residual_sum。 |
| L88–90 · 生成并汇总 `income_residual` | 计算每条观测的 income × residual,再把乘积之和保存到 income_residual_sum。 |
| L91–92 · `display ... scalar(...)` | 按 %21.12g 打印两个标量。这里的 % 格式控制显示,不改变原数值。 |
| L93 · `drop income_residual` | 删除仅用于核对的乘积变量,保留其汇总标量。 |
| L94–97 · 保存 / 导出 | 把含 fitted、residual、residual_sq 的 98 条明细写入 households_fitted.dta 和 CSV。 |
> 含截距 OLS 的正规方程:Σêᵢ = 0,ΣXᵢêᵢ = 0(数值计算允许微小舍入误差)。
现有输出:Σêᵢ = -1.955414e-11,ΣXᵢêᵢ = -7.333438e-08,都在浮点精度下接近 0。以 C4 为例,income=800、consume=561,拟合消费约 605.969927,残差约 -44.969927 元,表示实际消费低于拟合值。
> **不是模型假设检验**
>
> 这两个零和性质由当前含截距 OLS 的求解方式保证,不能证明总体 E(u|X)=0,不能证明不存在遗漏变量,也不能证明同方差。它们适合检查计算过程是否一致。
本段 meanonly 的行为依据 [S5];残差定义及 OLS 结果依据 [S7]。上面两个等式可直接由残差平方和对截距和斜率求导得到。
<a id="group"></a>
## 09 · 按收入水平生成汇总表
**原文件第 99–109 行**
```stata
preserve
collapse (count) n=consume (mean) mean_consume=consume (sd) sd_consume=consume (min) min_consume=consume (max) max_consume=consume, by(income)
recast double mean_consume sd_consume
generate double fitted_mean = scalar(intercept_ols) + scalar(slope_ols)*income
generate double mean_residual = mean_consume - fitted_mean
save "income_group_summary.dta", replace
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "income_group_summary.csv", datafmt replace
list, noobs separator(0)
restore
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L99 · `preserve` | 保存逐户明细,以便汇总后恢复。 |
| L100 · `collapse ... , by(income)` | 按收入分组,把 98 条明细替换成 10 条分组记录;(count) 数非缺失消费条数,(mean) 求均值,(sd) 求组内样本标准差,(min)/(max) 求范围。 |
| L100 · `n=consume` 等写法 | 左侧是汇总后的新变量名,右侧是被汇总的原变量。例如 mean_consume=consume 生成消费均值。 |
| L101 · `recast double ...` | 保证均值、标准差的后续存储类型为 double;类型扩展不能恢复此前若已发生的精度损失。 |
| L102 · `fitted_mean` | 用此前保存的截距和斜率计算每个收入水平的拟合均值;没有重新对 10 个均值做回归。 |
| L103 · `mean_residual` | 组内实际均值减去线性拟合值。同组收入相同,所以也等于该组逐户残差的平均值。 |
| L104–108 · 保存、导出并 `list` | 将分组表保存为 DTA/CSV,并在日志中显示。 |
| L109 · `restore` | 恢复 98 条逐户明细;后面的散点图因此基于家庭记录。 |
| 收入 | n | 组内均值 | 组内标准差 | 拟合均值 | 平均残差 |
| --- | --- | --- | --- | --- | --- |
| 800 | 4 | 605.000 | 34.785 | 605.970 | -0.970 |
| 1100 | 6 | 825.000 | 121.698 | 824.610 | +0.390 |
| 1400 | 11 | 1045.000 | 116.309 | 1043.249 | +1.751 |
| 1700 | 13 | 1265.000 | 137.683 | 1261.889 | +3.111 |
| 2000 | 13 | 1485.000 | 135.988 | 1480.528 | +4.472 |
| 2300 | 13 | 1682.154 | 154.879 | 1699.168 | -17.014 |
| 2600 | 13 | 1925.000 | 165.183 | 1917.808 | +7.192 |
| 2900 | 10 | 2136.200 | 125.612 | 2136.447 | -0.247 |
| 3200 | 9 | 2350.333 | 203.426 | 2355.087 | -4.754 |
| 3500 | 6 | 2583.167 | 248.066 | 2573.727 | +9.440 |
OLS 不要求每个收入组的平均残差都等于 0;这里满足的是全样本正规方程。各组人数不同,不能用 10 个组均值的不加权回归替代当前的逐户回归。组内标准差也不是均值的标准误。collapse 命令依据 [S10]。
<a id="plot"></a>
## 10 · 绘制散点与拟合直线,导出四种图形
**原文件第 111–117 行**
```stata
set scheme s2color
graph set window fontface "SimSun"
twoway (scatter consume income, mcolor("26 71 111") msymbol(O) msize(small)) (lfit consume income, lcolor("174 44 51") lwidth(medthin)), title("不同可支配收入水平组家庭消费支出的条件分布图", size(medium) color(navy)) xtitle("每月可支配收入(元)", size(medsmall)) ytitle("每月消费支出(元)", size(medsmall)) xlabel(1000(500)3500) ylabel(500(500)3000, angle(90) grid glcolor("225 236 238")) xscale(range(500 4000)) yscale(range(0 3600)) legend(order(1 "consume" 2 "Fitted values") rows(1) region(lcolor(black) fcolor(white))) graphregion(color("234 242 243")) plotregion(color(white)) xsize(9) ysize(6.55) name(consumption_fit, replace)
graph save "consumption_fit.gph", replace
graph export "consumption_fit.png", width(1800) replace
graph export "consumption_fit.svg", replace
graph export "consumption_fit.pdf", replace
```
| 行号 / 选项 | 解释 |
| --- | --- |
| L111 · `set scheme s2color` | 使用 s2color 图形样式方案。 |
| L112 · `graph set window fontface "SimSun"` | 把图形窗口字体设为宋体,以支持这里的中文标题与坐标轴文字;实际显示取决于本机字体。 |
| L113 · `twoway (...) (...)` | 把两个图层叠加在同一坐标系:scatter 散点层与 lfit 线性拟合层。 |
| `scatter consume income` | 先写纵轴变量消费,后写横轴变量收入;每条家庭记录对应一个点,完全相同的点会重合。 |
| `mcolor("26 71 111") msymbol(O) msize(small)` | 散点为指定 RGB 深蓝色、圆形符号、小尺寸。 |
| `lfit consume income` | 使用当前数据拟合带截距的一元 OLS 直线并绘制,不是把相邻散点连起来,也不直接读取 fitted 列。当前数据和设定相同,因此与前面的 OLS 拟合线一致。 |
| `lcolor("174 44 51") lwidth(medthin)` | 拟合线为指定 RGB 红色、中偏细线宽。 |
| `title(...) xtitle(...) ytitle(...)` | 分别设置主标题、横轴标题、纵轴标题及文字大小;navy 指定标题颜色。 |
| `xlabel(1000(500)3500)` | 横轴刻度从 1000 到 3500,间隔 500;这是刻度控制,不是筛选收入。 |
| `ylabel(500(500)3000, angle(90) grid glcolor(...))` | 纵轴刻度从 500 到 3000、间隔 500;标签旋转 90 度,添加指定浅色网格。 |
| `xscale(range(500 4000)) yscale(range(0 3600))` | 设置坐标轴范围,留出边距;不删除范围之外的记录,也不改变回归样本。 |
| `legend(order(1 ... 2 ...) rows(1) region(...))` | 图例按散点、拟合线的顺序排列成一行,并设黑色边框、白色背景。 |
| `graphregion(...) plotregion(...)` | 整张图背景设为浅蓝色,内部绘图区设为白色。 |
| `xsize(9) ysize(6.55)` | 指定图形物理宽高(英寸),决定总体比例;PNG 像素宽度在导出时另设。 |
| `name(consumption_fit, replace)` | 给内存图形命名,并允许替换同名内存图;与磁盘文件覆盖是不同层次。 |
| L114 · `graph save ...gph` | 保存 Stata 原生图形,可在 Stata 中重新打开。 |
| L115 · `graph export ...png, width(1800)` | 导出宽 1800 像素的位图。 |
| L116–117 · SVG、PDF 导出 | 导出可缩放矢量图与适合打印的 PDF;这些 replace 允许覆盖对应磁盘文件。 |
![现有 Stata 输出:家庭消费散点与 OLS 拟合线(2026-10-23)](consumption_fit.png)
这幅图没有置信带,也没有预测区间。固定收入处的竖直点列是该收入组的经验消费分布,不意味着已知真实总体分布。原文件 L113 是一整条长命令,本页面仅自动折行显示,没有改变原始语句。[S11]
<a id="robust"></a>
## 11 · 计算异方差稳健标准误并导出
**原文件第 119–136 行**
```stata
* Heteroskedasticity-robust standard errors as a sensitivity check.
regress consume income, vce(robust)
estimates store Robust
estimates save "consumption_robust.ster", replace
tempname robustcoef
postfile `robustcoef' str12 term double estimate double se double t double p double ci_low double ci_high using "robust_coefficients.dta", replace
foreach term in income _cons {
scalar tval = _b[`term'] / _se[`term']
scalar crit = invttail(e(df_r), .025)
post `robustcoef' ("`term'") (_b[`term']) (_se[`term']) (tval) (2*ttail(e(df_r), abs(tval))) (_b[`term']-crit*_se[`term']) (_b[`term']+crit*_se[`term'])
}
postclose `robustcoef'
preserve
use "robust_coefficients.dta", clear
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "robust_coefficients.csv", datafmt replace
restore
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L119 · `* ...` | 行首星号表示注释,说明下一段是稳健标准误的敏感性检查,不执行运算。 |
| L120 · `regress ... , vce(robust)` | 同一模型仍以 OLS 求系数,用异方差稳健协方差估计计算标准误;本例是非聚类 HC1,带 N/(N−k)=98/96 修正。 |
| L121–122 · `estimates store/save` | 分别在内存和磁盘保存稳健版本,名字为 Robust、文件为 consumption_robust.ster。 |
| L123–130 · 生成稳健系数表 | 复用 L71–78 的流程,但此时 _se[] 来自稳健协方差矩阵;脚本仍用 e(df_r)=96 的 t 分布计算 p 值和区间。 |
| L131–136 · 导出并恢复数据 | 将 robust_coefficients.dta 转成 CSV,最后恢复逐户明细。 |
> HC1 协方差 = [N/(N−k)] × (X′X)⁻¹ [Σ xᵢxᵢ′êᵢ²] (X′X)⁻¹;这里 k=2,X 含截距列。
| 变量 | 系数(两次相同) | 常规 SE | HC1 SE | 稳健 p 值 | 稳健 95% 区间 |
| --- | --- | --- | --- | --- | --- |
| income | 0.728799 | 0.020302 | 0.021428 | 2.38147e-55 | [0.686265, 0.771332] |
| _cons | 22.930934 | 47.096629 | 43.447563 | 0.598867 | [-63.311797, 109.173665] |
> **稳健标准误不一定更大**
>
> 本例收入系数的标准误从约 0.020302 增到 0.021428,而截距的标准误从约 47.096629 降到 43.447563。稳健处理针对方差估计,不会消除异方差,也不会解决遗漏变量或内生性。
这里不是抗异常值的稳健回归算法(例如 rreg),也没有使用聚类标准误。稳健 t 推断仍依赖相应独立性、外生性和大样本近似条件。HC1 修正及当前 regress 行为依据 Stata 15 自带官方帮助 [S7]。
<a id="hettest"></a>
## 12 · 恢复 OLS 并进行异方差检验
**原文件第 137–140 行**
```stata
estimates restore OLS
estat hettest
scalar bp_chi2 = r(chi2)
scalar bp_p = r(p)
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L137 · `estimates restore OLS` | 恢复 L56 保存的常规 OLS 估计结果。L136 的 restore 只恢复了数据,没有把当前稳健估计切回常规 OLS。 |
| L138 · `estat hettest` | 执行 Breusch–Pagan / Cook–Weisberg 异方差检验。无变量列表时,默认以拟合值为方差解释变量;默认 normal 版本有独立正态扰动假定。 |
| L139–140 · 保存 `r(chi2)` 和 `r(p)` | 把检验统计量及 p 值保存为标量,供最终模型汇总表使用。这里是 r() 而不是 e()。 |
现有结果:χ²(1) = 8.700549,p = 0.003181142。在默认检验条件下,5% 水平拒绝同方差原假设,提示方差可能随拟合水平变化;不能据此断定是哪条记录有问题。
> **这一行 restore 有实际作用**
>
> 本机 Stata 15 的 hettest.ado 会拒绝直接接在 vce(robust) 的结果后运行。因此 estimates restore OLS 是必要的模型切换;即使两次系数相同,协方差类型也不同。
不要把这里默认 normal 版本的统计量一概写成辅助回归的 N×R²;Stata 的 iid 选项才使用相应的 N×R² 版本。本文件没有添加 iid、fstat 或 rhs。命令行为已核对本机 regress_postestimation.sthlp 及 hettest.ado。[S7]
<a id="model-export"></a>
## 13 · 将模型汇总指标写成一行 CSV
**原文件第 142–161 行**
```stata
preserve
clear
set obs 1
generate N = scalar(N_ols)
generate intercept = scalar(intercept_ols)
generate slope = scalar(slope_ols)
generate r_squared = scalar(rsquared_ols)
generate root_mse = scalar(rmse_ols)
generate rss = scalar(rss_ols)
generate F = scalar(F_ols)
generate df_resid = scalar(df_ols)
generate p_F = scalar(pF_ols)
generate residual_sum = scalar(residual_sum)
generate income_residual_sum = scalar(income_residual_sum)
generate bp_chi2 = scalar(bp_chi2)
generate bp_p = scalar(bp_p)
quietly ds, has(type numeric)
format `r(varlist)' %21.15g
export delimited using "model_statistics.csv", datafmt replace
restore
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L142–144 · `preserve`、`clear`、`set obs 1` | 先保存当前家庭数据,再清空数据并创建 1 条空观测。这里是 clear,不是开头的 clear all;前面保存的标量仍可使用。 |
| L145–148 · `generate N/intercept/slope/r_squared` | 把样本量、截距、斜率、R² 标量复制为单行数据表的列。 |
| L149–153 · `root_mse/rss/F/df_resid/p_F` | 复制误差尺度、残差平方和、常规 F、残差自由度及常规 F 检验 p 值。F 取自最初 OLS 标量,不是稳健回归的 F。 |
| L154–157 · 残差和与 BP 结果 | 复制两项正规方程核对值及异方差检验统计量和 p 值。 |
| L158–160 · 格式与导出 | 再次选择数值变量并设格式,将 13 列、1 行的结果导出为 model_statistics.csv。原代码没有为这个汇总表另存 DTA。 |
| L161 · `restore` | 恢复原来的 98 条家庭明细。 |
generate 没有逐一写 double,是因为 L4 已设置默认类型;例如 generate N = scalar(N_ols) 得到一个数值变量。这里 N 只是列名,区别于 Stata 系统量 _N。该段的用意是把分散的标量整理成通用软件易于读取的一行表。[S1]、[S3]、[S6]
<a id="finish"></a>
## 14 · 完成标记、关闭日志与退出
**原文件第 163–165 行**
```stata
display "WEEK7_ANALYSIS_COMPLETED"
log close
exit, clear
```
| 行号 / 命令 | 解释 |
| --- | --- |
| L163 · `display "WEEK7_ANALYSIS_COMPLETED"` | 打印完成标记,用于定位是否执行到末尾。检查日志时应查找实际输出的标记行,而不只搜索被回显的 display 命令。 |
| L164 · `log close` | 关闭显式打开的分析日志,让日志完整写入磁盘。 |
| L165 · `exit, clear` | 在 do 文件中结束当前执行并返回调用者;clear 允许在退出 Stata 的情境下放弃未保存的数据,不是保存命令。通过 /e 启动时,批处理入口会在 do 执行结束后退出应用;若要从 do 内明确强制关闭 Stata,使用 exit, STATA clear。 |
既有 week7_stata.log 中有独立的完成标记输出行和关闭日志记录,没有 r(数字); 错误结束行。为核对退出语境,另在给定 Stata 可执行文件中运行了不读取原数据的最小测试:子 do 执行 exit, clear 后,父 do 的下一条 display 仍能执行;这个测试没有重跑原分析脚本。有关调用情境参见本机官方 [P] exit 与 [R] exit 帮助。[S12]
<a id="files"></a>
## 文件索引:每个输出保存了什么
| 文件 | 内容 / 规模 | 产生位置 |
| --- | --- | --- |
| `source_totals_audit.dta` / `.csv` | 10 个收入组的明细计数与合计差额 | L15–34 |
| `households_clean.dta` / `.csv` | 98 条家庭记录,4 个基础变量 | L49–52 |
| `consumption_ols.ster` | 常规 OLS 估计结果 | L57 |
| `regression_coefficients.dta` / `.csv` | 收入与截距两行常规系数表 | L72–84 |
| `households_fitted.dta` / `.csv` | 98 条记录;4 个基础变量 + fitted、residual、residual_sq | L94–97 |
| `income_group_summary.dta` / `.csv` | 10 组;收入、人数、均值、标准差、最小值、最大值、拟合均值、平均残差 | L104–107 |
| `consumption_fit.gph` / `.png` / `.svg` / `.pdf` | 同一散点与拟合线图的四种格式 | L114–117 |
| `consumption_robust.ster` | 稳健协方差版本的估计结果 | L122 |
| `robust_coefficients.dta` / `.csv` | 两行 HC1 系数推断表 | L124–136 |
| `model_statistics.csv` | 1 行 13 列的模型汇总,未另存 DTA | L160 |
| `week7_stata.log` | 脚本主动打开的纯文本日志 | L8、L164 |
| `week7_analysis.log` | 既有 /e 批处理执行产生的外层日志,不由脚本 log using 命名 | 运行方式产生 |
本次新建的 Week7_Stata代码详解.md 和 Week7_Stata代码详解.html 是代码解释文档;原 do 文件不会自动生成它们。原有“Week7_回归分析报告”是结果报告,与本次逐段代码说明分别保留。
<a id="checkpoints"></a>
## 理解与复现时最容易出错的地方
| 容易误读的地方 | 正确读法 |
| --- | --- |
| 把第 17 行合计当作额外家庭 | 合计只做审计;keep in 2/14 后它已被移除。 |
| 认为检出合计差额就已更正数据 | 程序只记录差异,回归始终使用可见明细。 |
| 认为所有 assert 都是统计假设检验 | 这里的 assert 是代码条件检查,失败时终止;它们没有检验显著性水平。 |
| 将空白当 0 或自动删掉重复金额 | 只删除 missing(consume);0 值和相同金额并不会因此被删掉。 |
| 把 excel_row 直接理解成家庭编号 | 它追踪原表行,长表中会跨收入重复;household_id 也是新建顺序号。 |
| 以为 format 可提高存储精度 | double/recast 管存储,format 管显示,datafmt 决定 CSV 是否沿用该显示格式。 |
| 在系数导出循环中插入另一次 regress | 会替换 e()、_b[]、_se[],使表格混入另一个模型结果。 |
| 将残差零和当作外生性证据 | 零和是含截距 OLS 的代数性质,不检验 E(u\|X)=0。 |
| 将图中的线当作各组均值连线 | lfit 是全样本一元 OLS 直线,组均值不必在线上。 |
| 以为 vce(robust) 消除了异方差 | 它调整推断所用协方差;当前系数、拟合值、残差不变。 |
| 把默认 hettest 当作无任何分布假定的检验 | 这里默认 normal 版本有相应假定,脚本没有请求 iid 版本。 |
| 把高 R² 解读成收入的因果作用 | 这是本组明细的单变量线性拟合,需要额外识别假设才能作因果解释。 |
<a id="run"></a>
## 如何在给定路径下运行和查阅
完整复现时,在 Stata 命令窗口运行下面的 do 命令。现有脚本包含 clear all、多处 replace 和末尾 exit;它面向完整批处理流程,重跑会覆盖同名分析输出。本文只说明原代码的这些实际行为。
```stata
do "C:/Users/cattom/Documents/econometrics/Week7/week7_analysis.do"
```
也可在 Windows PowerShell 启动独立 Stata 进程并等待完成。下面的可执行文件路径以 .exe 结束,不带末尾目录分隔符。此命令是供复现的示例,本次未执行。
```powershell
Start-Process -FilePath "D:\Stata15\StataSE-64.exe" `
-ArgumentList '/e do "C:\Users\cattom\Documents\econometrics\Week7\week7_analysis.do"' `
-WorkingDirectory "C:\Users\cattom\Documents\econometrics\Week7" `
-WindowStyle Hidden -Wait
```
读已有结果时,CSV 可以直接打开;DTA 用 use 读取,STER 用 estimates use 读取。例如:
```stata
use "C:/Users/cattom/Documents/econometrics/Week7/households_fitted.dta", clear
estimates use "C:/Users/cattom/Documents/econometrics/Week7/consumption_ols.ster"
estimates replay
```
若要交互式逐段学习,可在另存的练习副本中分段运行,并把宏定义、循环以及对应的 preserve/restore 放在同一执行块内。不要直接把带 L 行号的展示文本当作可执行命令;网页代码块的“复制代码”按钮复制的是无行号的原代码。
<a id="sources"></a>
## 核对依据与文档范围
以下来源均是本机可直接核查的文件,命令说明优先对应 Stata 15 安装附带的官方帮助。本文不依赖在线脚本、字体或公式服务。网页图形内嵌,可离线单独阅读;Markdown 使用同目录 consumption_fit.png 的相对图片引用。
| 编号 | 本机官方文档(相对 D:/Stata15/ado/base/) | 支撑内容 |
| --- | --- | --- |
| [S1] | `v/version.sthlp`、`c/clear.sthlp`、`g/generate.sthlp`、`m/macro.sthlp` | 版本控制、清空数据、新变量类型、局部宏等。 |
| [S2] | `p/post.sthlp` | postfile/post/postclose、文件句柄、逐行写入表。 |
| [S3] | `f/format.sthlp`、`i/import_delimited.sthlp` | 显示格式、datafmt 与 CSV 导出。 |
| [S4] | `i/import_excel.sthlp`、`r/reshape.sthlp` | Excel 导入与宽长表转换。 |
| [S5] | `s/summarize.sthlp` | summarize、meanonly、r(N)/r(sum) 等返回值。 |
| [S6] | `p/preserve.sthlp` | 保存并恢复当前数据集。 |
| [S7] | `r/regress.sthlp`、`r/regress_postestimation.sthlp`、`h/hettest.ado` | OLS、HC1、predict、e()、BP/CW 默认设定及 robust 后的限制。 |
| [S8] | `e/estimates_store.sthlp`、`e/estimates_save.sthlp` | 内存与磁盘估计结果存储。 |
| [S9] | `f/f_ttail.ihlp`、`f/f_invttail.ihlp`、`f/f_ftail.ihlp` | ttail、invttail、Ftail 的官方函数说明(由同名 .sthlp 引用)。 |
| [S10] | `c/collapse.sthlp` | 组内 count/mean/sd/min/max 汇总。 |
| [S11] | `t/twoway_lfit.sthlp`、`s/scatter.sthlp` 及相应 graph 帮助 | 散点、线性拟合层及图形选项。 |
| [S12] | `e/exit.sthlp`、`e/exit_program.sthlp` | 退出应用与退出程序/do 文件的语境。 |
数值依据:原始 Excel、week7_stata.log、households_clean.csv、households_fitted.csv、source_totals_audit.csv、regression_coefficients.csv、robust_coefficients.csv、income_group_summary.csv、model_statistics.csv。来源单元格、计数、系数、两种标准误、p 值与区间、拟合值、残差、R²、Root MSE、F、分组统计及默认 BP/CW 检验均已独立核对。
代码来源 SHA-256:`6b47d5ef559a1a68d866d73af54001ac84fb8d599debd4074b7858d22c673ffd`。行号与本文解释对应 2026 年 10 月 26 日读取的这份 165 行脚本;以后增删代码会改变行号。