原文链接：http://tecdat.cn/?p=20613

原文出处：拓端数据部落公众号

洛伦兹曲线来源于经济学，用于描述社会收入不均衡的现象。将收入降序排列，分别计算收入和人口的累积比例。
本文，我们研究收入和不平等。我们从一些模拟数据开始

1. 


2. > (income=sort(income))

3. [1]  19246  23764  53237  61696 218835

为什么说这个样本中存在不平等？如果我们看一下最贫穷者拥有的财富，最贫穷的人（五分之一）拥有5％的财富；倒数五分之二拥有11％，依此类推

1. > income[1]/sum(income)

2. [1] 0.0510

3. > sum(income[1:2])/sum(income)

4. [1] 0.1140

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

1. 


2. > plot(Lorenz(income))

3. > points(c(0:5)/5,c(0,cumsum(income)/sum(income))

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

1. > summary(income)

2. Min. 1st Qu.  Median    Mean 3rd Qu.    Max.

3. 2191   23830   42750   77010   87430 2003000

4. > hist(log(income),

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

1. > gini=function(x){

2. 


3. + mu=mean(x)

4. + g=2/(n*(n-1)*mu)*sum((1:n)*sort(x))-(n+1)/(n-1)

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

1. 


2. > G=boot(income,gini,1000)

3. > hist(G,col="light blue",border="white"

红色部分是90％置信区间，

1. 


2. 5%       95%

3. 0.4954235 0.5743917

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

1. > segments(quantile(G,.05),1,quantile(G,.95),1,

2. 


3. > lines(u,dnorm(u,mean(G),sd(G)),

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

1. 


2. > plot(x,y)

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

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

1. 


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

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

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

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

1. 


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

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

1. 


2. shape           rate

3. 1.0812757769   0.0140404379

4. (0.0604530180) (0.0009868055)

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

1. 


2. shape           rate

3. 1.0812757769   0.0014040438

4. (0.0473722529) (0.0000544185)

5. meanlog       sdlog

6. 6.11747519   1.01091329

7. (0.04520942) (0.03196789)

我们可以可视化密度

1. > hist(income,breaks=seq(0,2005000,by=5000),

2. + col=rgb(0,0,1,.5),border="white",

3. + fit_g$estimate[2])/1e2

4. + fit_ln$estimate[2])/1e2

5. > lines(u,v_g,col="red",lwd=2)

6. > lines(u,v_ln,col=rgb(1,0,0,.4),lwd=2)

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

1. > plot(x,y,type="s",col="black",xlim=c(0,250000))

2. 


3. + fit_g$estimate[2])

4. 


5. + fit_ln$estimate[2])

6. > lines(u,v_g,col="red",lwd=2)

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

对数据进行建模，

1. fit(ID=rep("Data",n),

2. 


3. 


4. 


5. Time difference of 2.101471 secs

6. for LNO fit across 1 distributions

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

1. > y2=N/sum(N)/diff(income_binned$low)

2. 


3. + fit_LN$parameters[2])

4. > plot(u,v,col="blue",type="l",lwd=2)

5. > for(i in 1:(n-1)) rect(income_binned$low[i],0,

6. + income_binned$high[i],y2[i],col=rgb(1,0,0,.2),

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

1. > v <- plnorm(u,fit_LN$parameters[1],

2. + fit_LN$parameters[2])

4. > for(i in 1:(n-1)) rect(income_binned$low[i],0,

5. 


6. 


7. > for(i in 1:(n-1)) rect(income_binned$low[i],

8. + y1[i],income_binned$high[i],c(0,y1)[i],

9. 


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

也可以拟合广义beta分布

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

为了获得最佳模型，查看

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

结果很好，接下来看下真实数据：

1. fit(ID=rep("US",n),

2. 


3. + distribution=LNO, distName="LNO"

4. Time difference of 0.1855791 secs

5. for LNO fit across 1 distributions

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

1. > v=dlnorm(u,fit_LN$parameters[1],

2. 


3. > plot(u,v,col="blue",type="l",lwd=2)

4. > for(i in 1:(n-1)) rect(data$low[i],

5. 


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

1. >

2. GB_family(ID=rep("US",n),

3. 


4. + ID_name="Country")

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

1. gini      aic      bic

2. 1  4.413431 825368.5 825407.3

3. 2  4.395080 825598.8 825627.9

4. 3  4.451881 825495.7 825524.8

5. 4  4.480850 825881.7 825910.8

6. 5  4.417276 825323.6 825352.7

7. 6  4.922122 832408.2 832427.6

8. 7  4.341036 827065.2 827084.6

9. 8  4.318667 826112.8 826132.2

10. 9        NA 831054.2 831073.6

11. 10       NA       NA       NA

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

最受欢迎的见解

1.R语言泊松Poisson回归模型分析案例

2.R语言进行数值模拟：模拟泊松回归模型

3.r语言泊松回归分析

4.R语言对布丰投针（蒲丰投针）实验进行模拟和动态可视化

5.用R语言模拟混合制排队随机服务排队系统

6.GARCH（1,1），MA以及历史模拟法的VaR比较

7.R语言做复杂金融产品的几何布朗运动的模拟

8.R语言进行数值模拟：模拟泊松回归模型

9.R语言对巨灾风险下的再保险合同定价研究案例：广义线性模型和帕累托分布Pareto distributions