R语言建模收入不平等:分布函数拟合及洛伦兹曲线(Lorenz curve)

简介: R语言建模收入不平等:分布函数拟合及洛伦兹曲线(Lorenz curve)

洛伦兹曲线来源于经济学,用于描述社会收入不均衡的现象。将收入降序排列,分别计算收入和人口的累积比例。

本文,我们研究收入和不平等。我们从一些模拟数据开始

> (income=sort(income))
[1]  19246  23764  53237  61696 218835

为什么说这个样本中存在不平等?如果我们看一下最贫穷者拥有的财富,最贫穷的人(五分之一)拥有5%的财富;倒数五分之二拥有11%,依此类推

> income[1]/sum(income)
[1] 0.0510
> sum(income[1:2])/sum(income)
[1] 0.1140


如果我们绘制这些值,就会得到 洛伦兹曲线

> plot(Lorenz(income))
> points(c(0:5)/5,c(0,cumsum(income)/sum(income))


现在,如果我们得到500个观测值。直方图是可视化这些数据分布的方法

> summary(income)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   2191   23830   42750   77010   87430 2003000 
> hist(log(income),


在这里,我们使用直方图将样本可视化。但不是收入,而是收入的对数(由于某些离群值,我们无法在直方图上可视化)。现在,可以计算 基尼系数 以获得有关不平等的一些信息

> gini=function(x){
 
+ mu=mean(x)
+ g=2/(n*(n-1)*mu)*sum((1:n)*sort(x))-(n+1)/(n-1)


实际上,没有任何置信区间的系数可能毫无意义。计算置信区间,我们使用boot方法

> G=boot(income,gini,1000)
> hist(G,col="light blue",border="white"


红色部分是90%置信区间,

5%       95% 
0.4954235 0.5743917


还包括了一条具有高斯分布的蓝线,

> segments(quantile(G,.05),1,quantile(G,.95),1,
 
> lines(u,dnorm(u,mean(G),sd(G)),


另一个流行的方法是帕累托图(Pareto plot),我们在其中绘制了累积生存函数的对数与收入的对数,

> plot(x,y)


如果点在一条直线上,则意味着可以使用帕累托分布来建模收入。

前面我们已经看到了如何获得洛伦兹曲线。实际上,也可以针对某些参数分布(例如,一些对数正态分布)获得Lorenz曲线,

> lines(Lc.lognorm,param=1.5,col="red")
> lines(Lc.lognorm,param=1.2,col="red")
> lines(Lc.lognorm,param=.8,col="red")


在这里, 对数正态分布是一个很好的选择。帕累托分布也许不是:

> lines(Lc.pareto,param=1.2,col="red")


实际上,可以拟合一些参数分布。

shape           rate    
  1.0812757769   0.0140404379 
 (0.0604530180) (0.0009868055)


现在,考虑两种分布,伽马分布和对数正态分布(适用于极大似然方法)

shape           rate    
  1.0812757769   0.0014040438 
 (0.0473722529) (0.0000544185)
     meanlog       sdlog   
  6.11747519   1.01091329 
 (0.04520942) (0.03196789)


我们可以可视化密度

> hist(income,breaks=seq(0,2005000,by=5000),
+ col=rgb(0,0,1,.5),border="white",
  + fit_g$estimate[2])/1e2
 + fit_ln$estimate[2])/1e2
> lines(u,v_g,col="red",lwd=2)
> lines(u,v_ln,col=rgb(1,0,0,.4),lwd=2)


在这里,对数正态似乎是一个不错的选择。我们还可以绘制累积分布函数

`> plot(x,y,type="s",col="black",xlim=c(0,250000))
+ fit_g$estimate[2])
+ fit_ln$estimate[2])
> lines(u,v_g,col="red",lwd=2)`


现在,考虑一些更现实的情况,在这种情况下,我们没有来自调查的样本,但对数据进行了合并,

对数据进行建模,

fit(ID=rep("Data",n), 
Time difference of 2.101471 secs
for LNO fit across 1 distributions


我们可以拟合对数正态分布(有关该方法的更多详细信息,请参见 从合并收入估算不平等 的方法)

`> y2=N/sum(N)/diff(income_binned$low) 
+ fit_LN$parameters[2])
> plot(u,v,col="blue",type="l",lwd=2)
> for(i in 1:(n-1)) rect(income_binned$low[i],0,
+ income_binned$high[i],y2[i],col=rgb(1,0,0,.2),`


在此,在直方图上(由于已对数据进行分箱,因此很自然地绘制直方图),我们可以看到拟合的对数正态分布很好。

> v <- plnorm(u,fit_LN$parameters[1],
+ fit_LN$parameters[2])
> for(i in 1:(n-1)) rect(income_binned$low[i],0,
> for(i in 1:(n-1)) rect(income_binned$low[i],
+ y1[i],income_binned$high[i],c(0,y1)[i],


对于累积分布函数,我考虑了最坏的情况(每个人都处于较低的收入中)和最好的情况(每个人都具有最高可能的收入)。

也可以拟合广义beta分布

GB_family(ID=rep("Fake Data",n),


为了获得最佳模型,查看

> fits[,c("gini","aic","bic")]


结果很好,接下来看下真实数据:

fit(ID=rep("US",n), 
+ distribution=LNO, distName="LNO"
Time difference of 0.1855791 secs
for LNO fit across 1 distributions


同样,我尝试拟合对数正态分布

> v=dlnorm(u,fit_LN$parameters[1],
> plot(u,v,col="blue",type="l",lwd=2)
> for(i in 1:(n-1)) rect(data$low[i],


但是在这里,拟合度很差。同样,我们可以估算广义beta分布

> 
GB_family(ID=rep("US",n), 
+ ID_name="Country")


可以得到基尼指数,  AIC 和BIC

gini      aic      bic
1  4.413431 825368.5 825407.3
2  4.395080 825598.8 825627.9
3  4.451881 825495.7 825524.8
4  4.480850 825881.7 825910.8
5  4.417276 825323.6 825352.7
6  4.922122 832408.2 832427.6
7  4.341036 827065.2 827084.6
8  4.318667 826112.8 826132.2
9        NA 831054.2 831073.6
10       NA       NA       NA


看到最好的分布似乎是 广义伽玛分布。


相关文章
|
27天前
|
数据挖掘 C语言 C++
R语言是一种强大的统计分析工具,提供了丰富的函数和包用于时间序列分析。
【10月更文挑战第21天】时间序列分析是一种重要的数据分析方法,广泛应用于经济学、金融学、气象学、生态学等领域。R语言是一种强大的统计分析工具,提供了丰富的函数和包用于时间序列分析。本文将介绍使用R语言进行时间序列分析的基本概念、方法和实例,帮助读者掌握R语言在时间序列分析中的应用。
42 3
|
6月前
【R语言实战】——带有高斯新息的金融时序的GARCH模型拟合预测及VAR/ES风险度量
【R语言实战】——带有高斯新息的金融时序的GARCH模型拟合预测及VAR/ES风险度量
|
6月前
【R语言实战】——带有新息为标准学生t分布的金融时序的GARCH模型拟合预测
【R语言实战】——带有新息为标准学生t分布的金融时序的GARCH模型拟合预测
|
6月前
|
数据可视化
【R语言实战】——金融时序分布拟合
【R语言实战】——金融时序分布拟合
|
6月前
|
机器学习/深度学习 数据可视化
R语言逻辑回归logistic模型ROC曲线可视化分析2例:麻醉剂用量影响、汽车购买行为2
R语言逻辑回归logistic模型ROC曲线可视化分析2例:麻醉剂用量影响、汽车购买行为
|
2月前
|
机器学习/深度学习
R语言模型评估:深入理解混淆矩阵与ROC曲线
【9月更文挑战第2天】混淆矩阵和ROC曲线是评估分类模型性能的两种重要工具。混淆矩阵提供了模型在不同类别上的详细表现,而ROC曲线则通过综合考虑真正率和假正率来全面评估模型的分类能力。在R语言中,利用`caret`和`pROC`等包可以方便地实现这两种评估方法,从而帮助我们更好地理解和选择最适合当前任务的模型。
|
3月前
|
数据采集 存储 数据可视化
R语言时间序列分析:处理与建模时间序列数据的深度探索
【8月更文挑战第31天】R语言作为一款功能强大的数据分析工具,为处理时间序列数据提供了丰富的函数和包。从数据读取、预处理、建模到可视化,R语言都提供了灵活且强大的解决方案。然而,时间序列数据的处理和分析是一个复杂的过程,需要结合具体的应用场景和需求来选择合适的方法和模型。希望本文能为读者在R语言中进行时间序列分析提供一些有益的参考和启示。
|
3月前
|
数据可视化 数据挖掘 数据处理
R语言函数与自定义函数:提高代码的复用性
【8月更文挑战第27天】 自定义函数是R语言编程中不可或缺的一部分,它们通过封装复杂的逻辑和提供灵活的参数化设计,极大地提高了代码的复用性和可维护性。通过掌握自定义函数的基本语法和高级技巧,我们可以编写出更加高效、可读的R语言代码,从而更好地应对复杂的数据分析和统计建模任务。
|
6月前
|
数据可视化
【R语言实战】——金融时序ARIMA建模
【R语言实战】——金融时序ARIMA建模
|
6月前
|
机器学习/深度学习 数据可视化
R语言Stan贝叶斯回归置信区间后验分布可视化模型检验|附数据代码
R语言Stan贝叶斯回归置信区间后验分布可视化模型检验|附数据代码