在R中计算基于分位数的分数

Question

我有一个数据框，包含2006年到2010年的4个行业部门、150家公司和这些公司的净收入。总共有750条观察记录，每年每家公司都有一条。我想根据五分位数为每个行业年度内的公司收入评分。所以，在每个行业年度内，净收入在前20%的公司得到5分，接下来20%得到4分，以此类推。最后20%得到1分。

样本数据库如下：

Year Industry Firm Income
2006 Chemicals ABC 334.50
2007 Chemicals ABC 388.98
.
.
2006 Pharma XYZ 91.45
.
.

我该如何在R中完成这个任务？我尝试过使用aggregate和tapply函数，以及quantile函数，但是我无法得出应该使用的逻辑。请帮忙。

我曾试图只将最低的20%分配为1分，但返回了错误。

db10$score <- ifelse(db10$income < aggregate(income~Year+industry,db10,quantile,c(0.2)),1,0)

- Piyush Shah

---

你能提供一个数据集的例子吗？ - desc

发布了数据集的典型样本。 - Piyush Shah

Answer 1

请尝试以下方法：

首先，我将创建一个示例来测试下面的函数:

y = c(rep(2001,15),rep(2002,15),rep(2003,15))
ind = c("A","B","C","D","E","G","H","I","J","K","L","M","N","O","P")
val = runif(45,10,100)
df = data.frame(y,ind,val)

head(df,20)

      y ind      val
1  2001   A 63.32011
2  2001   B 85.67976
3  2001   C 86.77527
4  2001   D 32.18319
5  2001   E 49.86626
6  2001   G 57.73214
7  2001   H 18.08216
8  2001   I 22.31012
9  2001   J 44.11174
10 2001   K 54.76902
11 2001   L 41.82495
12 2001   M 64.84514
13 2001   N 59.16529
14 2001   O 61.28870
15 2001   P 84.76561
16 2002   A 83.68185
17 2002   B 45.01354
18 2002   C 62.22964
19 2002   D 98.41717
20 2002   E 19.91548

有3个年份，行业从A到P。数据框按年和行业排序。

下面的函数接受一个年份值y并计算所有df$val的分位数类别，其中年份df$y是y

quintile = function(y) {
    x = df$val[df$y == y]
    qn = quantile(x, probs = (0:5)/5)
    result = as.numeric(cut(x, qn, include.lowest = T))
}

唯一剩下的事情就是将这个函数应用于唯一的年份值。

df$qn = unlist(lapply(unique(df$y), quintile))

结果：

> head(df,20)
      y ind      val qn
1  2001   A 63.32011  4
2  2001   B 85.67976  5
3  2001   C 86.77527  5
4  2001   D 32.18319  1
5  2001   E 49.86626  2
6  2001   G 57.73214  3
7  2001   H 18.08216  1
8  2001   I 22.31012  1
9  2001   J 44.11174  2
10 2001   K 54.76902  3
11 2001   L 41.82495  2
12 2001   M 64.84514  4
13 2001   N 59.16529  3
14 2001   O 61.28870  4
15 2001   P 84.76561  5
16 2002   A 83.68185  4
17 2002   B 45.01354  1
18 2002   C 62.22964  3
19 2002   D 98.41717  5
20 2002   E 19.91548  1

也许有更简单的方法来实现这个...

按两列分组

如果你想根据两列分组：y和grp来计算五分位数。

y = c(rep(2001,15),rep(2002,15),rep(2003,15))
grp = c("G1","G1","G1","G1","G1","G2","G2","G2","G2","G2","G3","G3","G3","G3","G3")
ind = c("A","B","C","D","E","G","H","I","J","K","L","M","N","O","P")
val = round(runif(45,10,100))
df = data.frame(y,grp,ind,val)

> head(df,20)
      y grp ind val
1  2001  G1   A  40
2  2001  G1   B  33
3  2001  G1   C  65
4  2001  G1   D  99
5  2001  G1   E  18
6  2001  G2   G  36
7  2001  G2   H  15
8  2001  G2   I  17
9  2001  G2   J  42
10 2001  G2   K  67
11 2001  G3   L  60
12 2001  G3   M  34
13 2001  G3   N  61
14 2001  G3   O  76
15 2001  G3   P  15
16 2002  G1   A  18
17 2002  G1   B  15
18 2002  G1   C  44
19 2002  G1   D  79
20 2002  G1   E  22

然后使用：

quintile = function(z) {
    x = df$val[df$y == z[1] & df$grp == z[2]]
    qn = quantile(x, probs = (0:5)/5)
    result = as.numeric(cut(x, qn, include.lowest = T))
}


df$qn = as.vector(apply(unique(df[,c("y","grp")]),1, quintile))

结果：

> head(df,20)
      y grp ind val qn
1  2001  G1   A  40  3
2  2001  G1   B  33  2
3  2001  G1   C  65  4
4  2001  G1   D  99  5
5  2001  G1   E  18  1
6  2001  G2   G  36  3
7  2001  G2   H  15  1
8  2001  G2   I  17  2
9  2001  G2   J  42  4
10 2001  G2   K  67  5
11 2001  G3   L  60  3
12 2001  G3   M  34  2
13 2001  G3   N  61  4
14 2001  G3   O  76  5
15 2001  G3   P  15  1
16 2002  G1   A  18  2
17 2002  G1   B  15  1
18 2002  G1   C  44  4
19 2002  G1   D  79  5
20 2002  G1   E  22  3

在这个例子中，y代表年份，grp代表行业组别，ind代表公司，val代表收入。

注意在apply中的c("y","grp")的顺序和quintile函数内部的列名顺序。你需要将它们替换为你想要的列名。

请注意，如果你的分组很小（在这个例子中每组只有5家公司），则五分位数可能不是唯一的，会出现错误提示。 使用问题中的列名

quintile = function(z) {
    x = df$Income[df$Year == z[1] & df$Industry == z[2]]
    qn = quantile(x, probs = (0:5)/5)
    result = as.numeric(cut(x, qn, include.lowest = T))
}


df$qn = as.vector(apply(unique(df[,c("Year","Industry")]),1, quintile))

在应用此操作之前，数据框 df 必须按照年份和行业进行排序。

- R. Schifini

---

@RSchifini 非常感谢，您提供的代码对您的数值有效。但是，很抱歉，我认为我没有表达清楚。我希望在每个行业年度内应用分位数得分。对于一年（比如2001年），行业A可能有许多公司（A1、A2、A3...）。评分需要在这些公司（A1、A2、A3...）内部工作，适用于每一年。 - Piyush Shah

@Rschifini我收到一个错误，上面写着Error in $<-.data.frame(tmp, "qn", value = list(1 = c(3, 5, 5, : replacement has 30 rows, data has 2398。我的数据有5年和6个组，所以才会出现30个数据的情况。我按照您的代码改变了变量名并运行了它。 - Piyush Shah

你的数据按年份和行业排序了吗？我会添加你的列名称代码。 - R. Schifini