该课程已经讲解完，感兴趣的小伙伴可以文末的阅读原文跳转到 RStata 平台观看录播～
在「经济评论」2022 年第 5 期的文章<数字普惠金融对农户收入不平等的影响研究———基于 CFPS 数据的实证分析>中使用了这两个指标：
Kakwani 指数：
Yitzhaki 指数：
在特定群组内，农户收入水平越高，受到的相对剥夺越低，即 Kakwani 指数越小，说明农户收入不平等程度越低。Yitzhaki 指数也类似，其值越小，农户收入不平等程度也越低。
下面我们一起来看下如何在 Stata 中计算这两个指标。
Kakwani 指数
首先我们读取数据，生成分组变量然后保存第一组的数据：
clear all
use test, clear
label drop _all
*- 按照年份和区县分组
egen group = group(year countyid)
keep if group == 1
list in 1/10
*> +---------------------------------------------+
*> | fid14 year countyid fincome1 group |
*> |---------------------------------------------|
*> 1. | 230033 2014 54 21900 1 |
*> 2. | 230034 2014 54 3350 1 |
*> 3. | 230035 2014 54 44000 1 |
*> 4. | 230037 2014 54 52600 1 |
*> 5. | 230038 2014 54 31000 1 |
*> |---------------------------------------------|
*> 6. | 230039 2014 54 20905 1 |
*> 7. | 230045 2014 54 22000 1 |
*> 8. | 230046 2014 54 34600 1 |
*> 9. | 230051 2014 54 3000 1 |
*> 10. | 230053 2014 54 72990 1 |
*> +---------------------------------------------+
Kakwani 指数的计算需要先构造相对剥夺矩阵 RD(yi, yj)，也就是下面代码中的 m1。例如 m1[1, 2] 表示的是个体 1 相对个体 2 被剥削的收入，如果个体 2 收入比个体 1 高，这个剥削就是 inc[2] - inc[1]，反之则是 0。
rowsum(m1) 就是计算个体 i 被相对其他个体被剥削的总和了。mean(inc) 是平均收入。rowsum(m1) / mean(inc) / n 就得到了 Kakwani 指数。
mata:
inc = st_data(., "fincome1")
n = rows(inc)
m1 = J(n, n, .)
for (i = 1; i <= n; i++) {
for (j = 1; j <= n; j++) {
if (inc[j, 1] - inc[i, 1] > 0) {
m1[i, j] = inc[j, 1] - inc[i, 1]
}
else {
m1[i, j] = 0
}
}
}
colsum(m1)
m_mean = rowsum(m1) / mean(inc) / n
stata("cap drop Kakwani_index")
st_addvar("double", "Kakwani_index")
st_store(., "Kakwani_index", m_mean)
end
最后几行是把计算得到的矩阵再保存到 dta 中：
*> +---------------------------------------------------------+
*> | fid14 year countyid fincome1 group Kakwani~x |
*> |---------------------------------------------------------|
*> 1. | 230033 2014 54 21900 1 .39599958 |
*> 2. | 230034 2014 54 3350 1 .89368602 |
*> 3. | 230035 2014 54 44000 1 .11828144 |
*> 4. | 230037 2014 54 52600 1 .05935769 |
*> 5. | 230038 2014 54 31000 1 .2452644 |
*> |---------------------------------------------------------|
*> 6. | 230039 2014 54 20905 1 .41872405 |
*> 7. | 230045 2014 54 22000 1 .3939441 |
*> 8. | 230046 2014 54 34600 1 .20415481 |
*> 9. | 230051 2014 54 3000 1 .90407761 |
*> 10. | 230053 2014 54 72990 1 0 |
*> +---------------------------------------------------------+
观察下收入和 Kakwani 指数的关系：
sc Kakwani_index fincome1
这个结果是符合预期的，也即取值在 0～1 之间、收入越高，Kakwani 指数越低。
Yitzhaki 指数
Yitzhaki 指数需要准备两个矩阵，一个是超过个体 i 的个体 j 的收入矩阵，另外一个是 0 1 矩阵，表示是否超过：
mata:
mata clear
inc = st_data(., "fincome1")
n = rows(inc)
m1 = J(n, n, .)
m2 = J(n, n, .)
for (i = 1; i <= n; i++) {
for (j = 1; j <= n; j++) {
if (inc[j, 1] - inc[i, 1] > 0) {
m1[i, j] = inc[j, 1]
m2[i, j] = 1
}
else {
m1[i, j] = 0
m2[i, j] = 0
}
}
}
m3 = rowsum(m2) / n
m_mean = rowsum(m1) :/ rowsum(m2) - inc
m4 = m3 :* m_mean
stata("cap drop Yitzhaki_index")
st_addvar("double", "Yitzhaki_index")
st_store(., "Yitzhaki_index", m4)
end
代码里面的 m1 是超过个体 i 的个体 j 的收入矩阵，计算行和就是计算超过 yi 的样本收入总和；m2 表示个体 i 的收入是否被个体 j 超过，所以 m2 的行和就是超过 yi 的样本数，两者对应元素相除以就得到了群组中收入超过个体 i 的样本的平均收入。再减去个体 i 的收入乘以百分比就得到 Yitzhaki 指数。
再检查下趋势：
sc Yitzhaki_index fincome1
同样符合预期。
循环计算所有的组
由于循环中不能出现 end，所以还是需要编写两个函数，然后在循环里调用：
*- 循环计算所有组的
clear all
use test, clear
label drop _all
*- 由于循环中不能出现 end，所以还是需要编写两个函数：
mata:
void Kakwani_index() {
inc = st_data(., "fincome1")
n = rows(inc)
m1 = J(n, n, .)
for (i = 1; i <= n; i++) {
for (j = 1; j <= n; j++) {
if (inc[j, 1] - inc[i, 1] > 0) {
m1[i, j] = inc[j, 1] - inc[i, 1]
}
else {
m1[i, j] = 0
}
}
}
colsum(m1)
m_mean = rowsum(m1) / mean(inc) / n
stata("cap drop Kakwani_index")
st_addvar("double", "Kakwani_index")
st_store(., "Kakwani_index", m_mean)
}
void Yitzhaki_index() {
inc = st_data(., "fincome1")
n = rows(inc)
m1 = J(n, n, .)
m2 = J(n, n, .)
for (i = 1; i <= n; i++) {
for (j = 1; j <= n; j++) {
if (inc[j, 1] - inc[i, 1] > 0) {
m1[i, j] = inc[j, 1]
m2[i, j] = 1
}
else {
m1[i, j] = 0
m2[i, j] = 0
}
}
}
m3 = rowsum(m2) / n
m_mean = rowsum(m1) :/ rowsum(m2) - inc
m4 = m3 :* m_mean
stata("cap drop Yitzhaki_index")
st_addvar("double", "Yitzhaki_index")
st_store(., "Yitzhaki_index", m4)
}
end
void 表示无返回值，我们这里并不需要返回值，只是希望调用函数的时候运行下函数里面的代码即可。
然后就可以循环计算各个组的了：
cap mkdir "res"
egen group = group(year countyid)
levelsof group, local(grouplist)
foreach g in `grouplist' {
di "`g'"
preserve
qui {
keep if group == `g'
mata: Kakwani_index()
mata: Yitzhaki_index()
save "res/`g'", replace
}
restore
}
再合并分组计算的结果：
local files: dir "res" files "*.dta"
di `"`files'"'
local firstfile: word 1 of `files'
use res/`firstfile', clear
drop in 1/`=_N'
foreach f in `files' {
append using "res/`f'"
}
save 计算结果, replace
*> +---------------------------------------------------------------------+
*> | fid14 year countyid fincome1 group Kakwani~x Yitzhak~x |
*> |---------------------------------------------------------------------|
*> 1. | 230034 2016 54 22100 10 .54683425 23226.154 |
*> 2. | 230035 2016 54 101500 10 .01883512 800 |
*> 3. | 230037 2016 54 111900 10 0 . |
*> 4. | 230038 2016 54 30000 10 .44305998 18818.462 |
*> 5. | 230045 2016 54 77000 10 .10757751 4569.2308 |
*> |---------------------------------------------------------------------|
*> 6. | 230046 2016 54 61000 10 .19450884 8261.5385 |
*> 7. | 230053 2016 54 8080 10 .80976529 34393.846 |
*> 8. | 230054 2016 54 43240 10 .3231672 13726.154 |
*> 9. | 230055 2016 54 24100 10 .51785714 21995.385 |
*> 10. | 230057 2016 54 19400 10 .59084323 25095.385 |
*> +---------------------------------------------------------------------+
最后在检查下两种指数和收入的关系：
tw sc Kakwani_index fincome1, ///
name(a, replace) nodraw ///
ti("Kakwani Index") xsc(log) ///
xla(100 1000 10000 100000 1000000)
tw sc Yitzhaki_index fincome1, ///
name(b, replace) nodraw ///
ti("Yitzhaki Index") xsc(log) ///
xla(100 1000 10000 100000 1000000)
gr combine a b, xsize(9.4) ysize(4)
都是符合预期的。
直播信息
该课程已经讲解完，感兴趣的小伙伴可以文末的阅读原文跳转到 RStata 平台观看录播～
直播地址：腾讯会议(需要报名 RStata 培训班参加) 讲义材料：需要购买 RStata 名师讲堂会员，详情可阅读：一起来学习 R 语言和 Stata 啦！学习过程中遇到的问题也可以随时提问！
更多关于 RStata 会员的更多信息可添加微信号 r_stata 咨询：
附件下载（点击文末的阅读原文即可跳转）：
https://rstata.duanshu.com/#/brief/course/713ae80a33b942fc9f83fba417bd014d
