clear all
set more off
set seed 20260601
global dir_in "/Users/ac/Desktop/使用 Stata 测算各区县鸟类丰度指数"
global dir_out "$dir_in/output"
use "$dir_out/panel_county_month_year_2015_2020.dta", clear
summarize year
local ymin = r(min)
local ymax = r(max)
display "面板: " _N " 县×月×年单元; 年份 " `ymin' " - " `ymax'
gen byte treat = inlist(city, "湖州市", "衢州市", "南昌市", "九江市",
"广州市", "贵阳市", "兰州市", "哈密市")
gen byte post = (year >= 2017)
count if treat == 1
display "处理组城市(县×月×年单元数): " r(N) " 对照组: " _N - r(N)
gen byte period = (year <= 2016)
preserve
collapse (mean) g = gstd_all [aw=n_checklist], by(treat period)
reshape wide g, i(treat) j(period)
gen diff = g0 - g1
summ diff if treat == 1, meanonly
local dt = r(mean)
summ diff if treat == 0, meanonly
local dn = r(mean)
local naive_did = `dt' - `dn'
display _n "朴素 DiD (丰度指数) = " %9.4f `naive_did'
display "论文表1 DID系数 β1 显著为正 → 若本指数方向一致，应 >0。"
local naive_did_s : display %9.4f `naive_did'
restore
preserve
collapse (mean) s = sp_richness [aw=n_checklist], by(treat period)
reshape wide s, i(treat) j(period)
gen diff = s0 - s1
summ diff if treat == 1, meanonly
local dt = r(mean)
summ diff if treat == 0, meanonly
local dn = r(mean)
local naive_did_sp = `dt' - `dn'
display _n "朴素 DiD (物种数 richness) = " %9.4f `naive_did_sp'
display "注意：本物种数为未做努力量校正原始计数，论文种类异质性一般也做 Liang 校正。"
local naive_did_sp_s : display %9.4f `naive_did_sp'
restore
preserve
collapse (mean) g = gstd_all [aw=n_checklist], by(treat year)
twoway
(connected g year if treat == 0, lcolor(gs12) lpattern(solid) mcolor(gs12) lwidth(*1.2))
(connected g year if treat == 1, lcolor("228 26 28") lpattern(dash) mcolor("228 26 28") lwidth(*1.8)),
legend(order(1 "非试验区" 2 "试验区") size(*0.85))
xline(2016.5, lpattern(dot) lcolor("228 26 28"))
xtitle("年份", size(small)) ytitle("相对丰度指数(标准化)均值", size(small))
title("图A 试验区 vs 非试验区 年度丰度轨迹", size(medium) color(black))
subtitle("数据处理 & 绘制：微信公众号 RStata", size(small))
note("红色虚线 = 2017 年政策实施；数据来源：RStata 数据中心（样本 2015–2020）", size(vsmall))
ylabel(, format(%9.2f)) xlabel(2015(1)2020, format(%9.0f))
graphr(margin(medium))
graph save `"$dir_out/_figA"', replace
restore
preserve
collapse (mean) water = gstd_waterbird (mean) land = gstd_landbird [aw=n_checklist], by(month)
twoway
(connected water month, lcolor("55 126 184") mcolor("55 126 184") lwidth(*1.5))
(connected land month, lcolor("228 26 28") mcolor("228 26 28") lwidth(*1.5)),
legend(order(1 "水鸟(多为候鸟)" 2 "陆鸟(多为留鸟)") size(*0.85))
xtitle("月份", size(small)) ytitle("相对丰度指数(标准化)均值", size(small))
title("图B 季节模式：水鸟冬季(11-2月)峰值", size(medium) color(black))
subtitle("数据处理 & 绘制：微信公众号 RStata", size(small))
note("与「候鸟冬季南下」生态事实一致；数据来源：RStata 数据中心（样本 2015–2020）", size(vsmall))
ylabel(, format(%9.2f)) xlabel(1(1)12, format(%9.0f))
graphr(margin(medium))
graph save `"$dir_out/_figB"', replace
restore
preserve
collapse (mean) gstd_waterfowl gstd_shorebird gstd_waterbird gstd_landbird [aw=n_checklist], by(year)
twoway
(connected gstd_landbird year, lcolor("228 26 28") mcolor("228 26 28") lwidth(*1.5))
(connected gstd_waterfowl year, lcolor("55 126 184") mcolor("55 126 184") lwidth(*1.5))
(connected gstd_shorebird year, lcolor("77 175 74") mcolor("77 175 74") lwidth(*1.5))
(connected gstd_waterbird year, lcolor("152 78 163") mcolor("152 78 163") lwidth(*1.5)),
legend(order(1 "陆鸟 landbirds" 2 "雁鸭类 waterfowl"
3 "鸻鹬鸥类 shorebirds" 4 "其他水鸟 waterbirds") size(*0.85))
xtitle("年份", size(small)) ytitle("相对丰度指数均值", size(small))
title("图C 四大分类群年度趋势", size(medium) color(black)) xlabel(2015(1)2020, format(%9.0f))
subtitle("数据处理 & 绘制：微信公众号 RStata", size(small))
note("数据来源：RStata 数据中心（样本 2015–2020）", size(vsmall))
ylabel(, format(%9.2f))
graphr(margin(medium))
graph save `"$dir_out/_figC"', replace
restore
gen region = "其他"
replace region = "东部" if prov == "北京市" | prov == "天津市" | prov == "河北省" | prov == "上海市" |
prov == "江苏省" | prov == "浙江省" | prov == "福建省" | prov == "山东省" |
prov == "广东省" | prov == "海南省"
replace region = "中部" if prov == "山西省" | prov == "安徽省" | prov == "江西省" | prov == "河南省" |
prov == "湖北省" | prov == "湖南省" | prov == "广西壮族自治区"
replace region = "西部" if prov == "内蒙古自治区" | prov == "重庆市" | prov == "四川省" | prov == "贵州省" |
prov == "云南省" | prov == "西藏自治区" | prov == "陕西省" | prov == "甘肃省" |
prov == "青海省" | prov == "宁夏回族自治区" | prov == "新疆维吾尔自治区"
replace region = "东北" if prov == "辽宁省" | prov == "吉林省" | prov == "黑龙江省"
preserve
collapse (mean) g = gstd_all [aw=n_checklist], by(region)
keep if inlist(region, "东部", "中部", "西部", "东北")
graph bar g, over(region, sort(g) label(angle(0) labsize(*0.8)))
ytitle("相对丰度指数均值", size(small))
title("图D 区域基线丰度(描述性)", size(medium) color(black))
subtitle("数据处理 & 绘制：微信公众号 RStata", size(small))
note("数据来源：RStata 数据中心（样本 2015–2020）", size(vsmall))
bar(1, color("55 126 184"))
ylabel(, format(%9.2f)) graphr(margin(medium))
graph save `"$dir_out/_figD"', replace
restore
graph combine `"$dir_out/_figA.gph"' `"$dir_out/_figB.gph"' `"$dir_out/_figC.gph"' `"$dir_out/_figD.gph"',
rows(2) cols(2)
title("与刘钊等 (2025) 论文趋势对照", size(medium) color(black))
subtitle("数据处理 & 绘制：微信公众号 RStata", size(small))
note("年度轨迹 / 季节模式 / 分类群趋势 / 区域基线；数据来源：RStata 数据中心「1980～2025 年观鸟记录、经纬度及其所处的省市区县数据」（样本 2015–2020）", size(vsmall))
graphr(margin(medium)) xsize(20) ysize(12)
graph export `"$dir_out/fig4_与刘钊论文对照_2015_2020.png"', width(4800) replace
erase `"$dir_out/_figA.gph"'
erase `"$dir_out/_figB.gph"'
erase `"$dir_out/_figC.gph"'
erase `"$dir_out/_figD.gph"'
display _n "已保存 $dir_out/fig4_与刘钊论文对照_2015_2020.png"
display _n "=== 描述性趋势汇总 ==="
preserve
collapse (mean) raw = gamma_hat [aw=n_checklist], by(year)
list year raw, sep(0) noobs
restore
display _n "分类群 2015 vs 2020 变化(标准化指数):"
preserve
keep if inlist(year, 2015, 2020)
collapse (mean) gstd_all gstd_waterbird gstd_landbird [aw=n_checklist], by(year)
list year gstd_all gstd_waterbird gstd_landbird, sep(0) noobs
restore
display _n "全部完成。"
|
评论