我心中的标杆

刚才我大表哥用QQ传给我一张今天拍的照片,一看是老家的,心里一阵激动。在外看了很多山山水水,到头来还是家乡美!

不过今天我主要不是想说家乡的风景多么清秀,而是想写一写一直想写但一直未写、自幼一起在这美丽土地长大的我的大表哥。

表哥上次答应给我小舅(不是他爸爸)物色一辆车,工作繁忙之余,雷厉风行地把车、手续什么的都弄好了,现在还一个人亲自千里迢迢从深圳开十几个小时送回了家,个中辛苦我很清楚,更重要的是,这是他又一次而且不知是第多少次不辞劳苦利索地完成他答应过别人的事情。

不夸张地说,大表哥是我心中一个做人的高标杆。我想用一句别人曾经写给我、我受之有愧而表哥当之无愧的话:无比靠谱,无比稳重,看着给人安心的感觉。言出必行是一个苛刻的道德标准,但就我所知,我难以在记忆中找出他食言的例子,哪怕有时候感觉他是随便说说的话,日后他都兑现了。

表哥在我多次的困难时刻,助我一臂之力、为我雪中送炭,否则我走不到今天。我真的欠他许多,而我一直只是默默地接受,从没有在他面前说出我的感激。但我心知肚明,我知道我该怎么做,怎么说已不重要。

 

统计那个显著性

最近在研究行业周期性的时候,用到了交叉相关系数(Cross-Correlation Function,CCF) 来考察某行业营收季度累计同比增速与GDP季度累计同比增速的滞后、同步或超前的相关性。CCF与自相关系数(Auto-Correlation Function,ACF)类似,只不过CCF是两个时间序列,而ACF是单个时间序列。

同事刚才问我:CCF(取了绝对值)高与CCF显著有什么区别?我说CCF不高未必不显著。上图中绝对值最大的CCF才0.4多一点,但它是显著的,因为它超出了CCF统计显著性的蓝色临界虚线,该虚线是由CCF统计显著性检验统计量(CCF的一种函数)的临界值计算出来的。

统计中的“显著”与现实生活中“明显”之类的说法有区别,但是我觉得它还是源于现实生活,某种意义上可以说是对现实生活中的感觉的一种理性量化。

上图中0.4几的相关系数似乎“不高”,但是它在我们假设的统计分布中却是“异常”的高。打一个不知恰不恰当的比方:在高中的时候,满分100分的英语考卷,我考了80分,也许很一般般,但是研究生入学考试中我的英语考了80分,那么就是很优异的成绩了,因为这是两种不同难度的英语考试,用统计的话说就是,这两种考试的成绩分布不同。

不知道我上面的理解对不对。统计中的很多概念我也有点晕,比如我很讨厌的第一类错误、第二类错误、势这些东东。

 

R的内存管理问题

昨天在公司的技术交流会上谈到了R的内存管理问题,有点意思,不过这个比较底层的东东我懂得还不多。

R用的是所谓的GC(Garbage Collection)的内存管理方式,它在程序后台自行管理内存,使编程人员摆脱一些计算机细节,而专注于程序的逻辑主体。目前使用这种方式的编程语言包括Java、C#、S(R属于S的变种)等。在这种方式下,应用程序中删除的变量不会立即释放内存空间,GC会自行计算最优的时机来进行统一的清理,举个简单的例子:

?Download eg1.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
> a = matrix(1:1e6, 1e3, 1e3)
> memory.size()
[1] 19.59
> b = a
> memory.size()
[1] 19.59
> b = b + 1
> memory.size()
[1] 27.28
> rm(a)
> memory.size()
[1] 27.28
> gc()
          used (Mb) gc trigger (Mb) max used (Mb)
Ncells  220720  5.9     467875 12.5   350000  9.4
Vcells 1256462  9.6    4628761 35.4  5223293 39.9
> memory.size()
[1] 19.1
> object.size(b)
8000112 bytes

我们在32位4G内存Win7操作系统上用R做千万次级别的显示循环(暂时没想到怎么非显式),该循环的功能就是反复对同一个变量赋值,但是循环到一定的时候操作系统会因内存耗尽而死机,原因可能是R的GC内存管理方式造成的,举个例子:

?Download eg2.r
1
2
3
4
5
6
7
8
9
10
> memory.size()
[1] 14.01
> system.time(for(i in 1:1e7) x = i ^ 2)
   user  system elapsed
   5.85    0.00    5.90
> memory.size()
[1] 49.89
> invisible(gc())
> memory.size()
[1] 8.76

对同一个变量反复赋值,内存消耗却不断变大,猜测是一些中间隐形的东东给占用了。一个非常笨的办法就是将这千万次循环拆分成很多个子循环,每个子循环完了之后再gc()一下。在不改变硬件设施并且时间不能慢死人的条件下,不知道哪位客官有没有更好的解决方案?

【附1】统计之都论坛对此问题的讨论:http://cos.name/cn/topic/106566
【附2】统计之都论坛对R快速运算的讨论::http://cos.name/cn/topic/106723

三八女生节怎么送礼物

晨读的时候在统计之都论坛看到了一个好玩的帖子,大概意思是说,有个哥们的班里要求明天男生给女生送礼物,男14人,女23人,每个男生要送两个女生礼物,个别男生送一个,但怎么分配成了一个问题……

男生、女生的学号如下:

?Download studentIDs.r
1
2
3
4
5
6
> boys = c(22, 29, 31, 32, 33, 35, 36, 37, 41, 42, 47, 52, 55, 56)
> boys
 [1] 22 29 31 32 33 35 36 37 41 42 47 52 55 56
> girls = setdiff(22:58, boys)
> girls
 [1] 23 24 25 26 27 28 30 34 38 39 40 43 44 45 46 48 49 50 51 53 54 57 58

某破解宋词密码的著名网友给出了一个简单的抽签办法:

?Download 38present.r
1
2
3
4
5
6
7
8
> set.seed(20120307)
> g.rearrange = sample(girls)
> b.rearrange = rep(sample(boys), length = length(girls))
> sort(paste(b.rearrange, "->", g.rearrange))
 [1] "22 -> 28" "29 -> 25" "29 -> 26" "31 -> 27" "31 -> 58" "32 -> 50" "33 -> 24"
 [8] "35 -> 30" "35 -> 46" "36 -> 45" "36 -> 51" "37 -> 23" "37 -> 57" "41 -> 54"
[15] "42 -> 38" "42 -> 40" "47 -> 39" "52 -> 49" "52 -> 53" "55 -> 44" "55 -> 48"
[22] "56 -> 34" "56 -> 43"

问题迎刃而解,统计学和R语言为女生节增光添彩,哈哈。

在微博上还看到一个极品故事:为迎接女生节,机电某班四十多个男生一人掏了一百块钱,给他们班唯一的一个女生买了部iPhone 4S。不过女生就应该这么宠!

另外,鄙人多年来对Women’s Day的中文翻译颇有微词,Women怎么就是“妇女”呢!我很生气,建议两会提议案更改名称。

行业分类

研究了一段时间的行业量化模型,行业分类是该模型的基石之一。行业分类的依据主要是公司收入或利润的来源比重。目前常用的行业分类有证监会行业分类、中证行业分类、GICS分类等。

证监会行业分类

顾名思义,是证监会发布的行业分类标准。分为13个门类、22个次类、91个大类、288个中类。由于该标准在国内的研究报告中使用较广,而且可以直接从我们公司的数据库中提取较详尽的数据,所以我们在模型中采纳了该分类标准。次类的数目比较适中,符合我们模型的设计要求。

中证行业分类

上海证券交易所与中证指数有限公司参照GICS分类编制的行业分类标准。该分类将上市公司分为10个一级行业和25个二级行业,没有再作更细致的划分。出于另外的考虑,我们的选股模型采用的是这个分类标准。今后我们会将行业模型和选股模型结合起来,这里就会涉及到证监会行业分类和中证行业分类的比较和统一问题。

GICS分类 

全球行业分类标准(Global Industry Classification Standard,GICS)是由标准普尔(S&P)与摩根斯坦利公司(MSCI)于1999年8月联手推出的行业分类系统。该标准为全球金融业提供了一个全面的、全球统一的经济板块和行业定义,已经在世界范围内得到广泛认可。它的意义在于不仅为创造易复制的、量体裁衣的投资组合提供了坚实基础,更使得对全球范围经济板块和行业的研究更具可比性。GICS为四级分类,包括10个行业类别,24个行业组,65个行业和146个子行业。国内的行业分类标准基本都借鉴了GISC标准,但比GICS标准更符合中国国情。

 

R小私集

今后会在这里持续收藏一些自编或自改的人民群众喜闻乐见的R小函数或语句,以备不时之需:

Fibonacci函数:

?Download fibonacci.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
fibonacci = function(n)
{
    ## Author: Ryan Rosario & Jason Fan
    ## Launch Date: 2010-07-27   
    ## Update Date: 2012-01-20
    ## Purpose: calculate Fibonacci sequence
 
    if(n < 1) stop("Input must be an integer >= 1")
    if (n == 1 | n == 2) 1
    else fibonacci(n - 1) + fibonacci(n - 2)
}
 
> sapply(1:20, fibonacci)
 [1]    1    1    2    3    5    8   13   21   34   55   89  144  233  377
[15]  610  987 1597 2584 4181 6765

读取Excel 2007文件的函数:

?Download readXLSX.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
readXLSX = function(file.path, sheet.name)
{
    ## Author: Jason Fan
    ## Launch Date: 2011-11-01   
    ## Update Date: 2012-01-10
    ## Purpose: read .xlsx format file
 
    library(RODBC)
    ch.xlsx = odbcConnectExcel2007(file.path)
    XLSX = sqlFetch(ch.xlsx, sheet.name, as.is = T)
    close(ch.xlsx)
    detach("package:RODBC")
    XLSX
}

读取数据库数据的函数:

?Download readDB.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
readDB = function(sql.code)
{
    ## Author: Edmund Lu & Jason Fan
    ## Launch Date: 2011-11-01   
    ## Update Date: 2012-01-10
    ## Purpose: read data from a database
 
    library(RODBC)
    ch.db = odbcConnect(dsn = '××××', uid = '××××',  pwd = '××××')
    DB = sqlQuery(ch.db, sql.code, as.is = T)
    close(ch.db)
    detach('package:RODBC')
    DB
}

反秩函数:

?Download revRank.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
revRank = function(x, na.last = 'keep', ties.method = 'average')
{
    ## Author: Jason Fan
    ## Launch Date: 2011-12-01   
    ## Update Date: 2012-02-20
    ## Purpose: opposite to rank() function
    ## Note: read rank() help file for arguments na.last and ties.method 
 
    rank(-x, na.last, ties.method)
}
 
>  set.seed(20120220)
>  x = sample((-2):7, replace = T)
>  x
 [1]  1  1 -1  5  6 -2  4  3 -1  0
>  rank(x)
 [1]  5.5  5.5  2.5  9.0 10.0  1.0  8.0  7.0  2.5  4.0
>  revRank(x)
 [1]  5.5  5.5  8.5  2.0  1.0 10.0  3.0  4.0  8.5  7.0
>  names(x) = letters[1:10]
>  x
 a  b  c  d  e  f  g  h  i  j 
 1  1 -1  5  6 -2  4  3 -1  0 
>  rank(x)
   a    b    c    d    e    f    g    h    i    j 
 5.5  5.5  2.5  9.0 10.0  1.0  8.0  7.0  2.5  4.0 
>  revRank(x)
   a    b    c    d    e    f    g    h    i    j 
 5.5  5.5  8.5  2.0  1.0 10.0  3.0  4.0  8.5  7.0

月跨度函数:

?Download monthSpan.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
monthSpan = function(from.date, to.date)
{
    ## Author: Jason Fan & Yishuo Deng
    ## Launch Date: 2011-09-15    
    ## Update Date: 2012-02-01
    ## Purpose: Count the number of months between 2 dates 
 
    date1 = paste(substr(from.date, 1, 7), "01", sep = "-")
    date2 = paste(substr(to.date, 1, 7), "01", sep = "-")
    number.of.months = length(seq(from = as.Date(date1), to = as.Date(date2), by = "month"))
    number.of.months
}
 
> monthSpan(from.date = "2010-03-25", to.date = "2010-06-28")
[1] 4

找重函数:

?Download findRep.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
findRep = function(x)
{
    ## Author: Jason Fan
    ## Launch Date: 2012-02-29   
    ## Update Date: 2012-02-29
    ## Purpose: find the repeated element(s) in a vector
 
    if(is.character(x)) names(table(x)[table(x) > 1]) 
    if(is.logical(x)) as.logical(names(table(x)[table(x) > 1])) else 
        as.numeric(names(table(x)[table(x) > 1]))
}
 
> set.seed(123)
> x1 = sample(10, replace = T)
> x2 = sample(letters, 10, replace = T)
> x3 = c(T, F, T)
> x1
 [1]  3  8  5  9 10  1  6  9  6  5
> x2
 [1] "y" "l" "r" "o" "c" "x" "g" "b" "i" "y"
> x3
[1]  TRUE FALSE  TRUE
> findRep(x1)
[1] 5 6 9
> findRep(x2)
[1] "y"
> findRep(x3)
[1] TRUE

集合运算函数:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
> (A <- c(sort(sample(1:20, 9)),NA))
 [1]  8 10 11 12 14 15 16 18 20 NA
> (B <- c(sort(sample(3:23, 7)),NA))
[1]  4  6  9 10 13 17 21 NA
> 
> # A∩B
> intersect(A, B)
[1] 10 NA
> 
> # A∪B 
> union(A, B)
 [1]  8 10 11 12 14 15 16 18 20 NA  4  6  9 13 17 21
> 
> # A-B & B-A
> setdiff(A, B)
[1]  8 11 12 14 15 16 18 20
> setdiff(B, A)
[1]  4  6  9 13 17 21
>  
> # A∪B = (A-B)∪(A∩B)∪(B-A)
> setequal(union(A, B), c(setdiff(A, B), intersect(A, B), setdiff(B, A)))
[1] TRUE

怎么使用tapply

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
> set.seed(20120507)
> FactorData = data.frame(stock_id = rep(1:4, 4), 
+                         report_date = rep(paste(rep(2008:2011, each = 4), "12-31", sep = "-")), 
+ net_profit = sample(1e8:2e8, 16), 
+ industry_character = rep(LETTERS[1:4], 4))
> FactorData
   stock_id report_date net_profit industry_character
1         1  2008-12-31  143053699                  A
2         2  2008-12-31  143445207                  B
3         3  2008-12-31  147623078                  C
4         4  2008-12-31  167865770                  D
5         1  2009-12-31  102574966                  A
6         2  2009-12-31  101809770                  B
7         3  2009-12-31  150244678                  C
8         4  2009-12-31  112340973                  D
9         1  2010-12-31  198702380                  A
10        2  2010-12-31  162393928                  B
11        3  2010-12-31  141121733                  C
12        4  2010-12-31  120825197                  D
13        1  2011-12-31  189642339                  A
14        2  2011-12-31  157312958                  B
15        3  2011-12-31  165180411                  C
16        4  2011-12-31  118997537                  D
> FactorMat = tapply(FactorData$net_profit, list(FactorData$report_date, FactorData$stock_id), function(x) x)
> FactorMat
                   1         2         3         4
2008-12-31 143053699 143445207 147623078 167865770
2009-12-31 102574966 101809770 150244678 112340973
2010-12-31 198702380 162393928 141121733 120825197
2011-12-31 189642339 157312958 165180411 118997537
> FactorSum = tapply(1:nrow(FactorData), list(FactorData$industry_character), function(x) sum(FactorData[x, "net_profit"]))
> FactorSum
        A         B         C         D 
633973384 564961863 604169900 520029477

第一个圣诞祝福

第一个圣诞祝福是我导师发过来的,这是真诚的祝福,我这么冷的心,意外之余还是不禁感动。

相处久了,就会磨出感情;离别久了,难免沉淀记念。

希望大家都好,那些喜乐幽怨已不再重要,哪怕今后物是人已非,记得我们曾经在一起。

 

投资建模小体会

玩了几个月的真数据,记录一些小经验,会持续更新:

  1. 教科书上面的数据不好玩。
  2. 时间序列分析学校里学得很没劲,但投资中很有味。
  3. 用历史数据建模,很容易犯将当时不知道的数据用到模型中的错误,所以得出来的结果可能会很好,用老板的话说就是:Too Good To Be True。
  4. 投资用的模型不追求严格的理论基础,有合理的解释即可,实用效果最重要。
  5. 投资模型实盘操作中犯了错就没有挽回的机会。
  6. 正确并且符合模型格式的数据的整理构造很耗时但很重要,否则就是垃圾进垃圾出,更严重的是丢失真金白银。
  7. 数据格式的规范统一会提高编程的效率。
  8. 程序变量的命名一定要直观,长一点不要紧,否则编写了几个月的程序,回过头来自个儿都看不懂。
  9. 建模的目的就是为了比市场早一点知道,早一点点就意味着胜出机会。
  10. 用最简单的模型获取最高的收益是我们不懈的追求。
  11. 发现模型结果不合逻辑但又不知错误在哪,是一件抓狂的事情,此所谓当局者迷。
  12. 写程序不能一次堆一个大蛋糕,而是要先做好一个一个的小蛋糕,然后拼成一个大蛋糕,也就是要模块化,并且每个模块不能太大。换句话说,模型其实是由很多功能相对独立的函数组装起来的。
  13. 模型不完善可以暂时容忍,但是错误必须杜绝,对模型的信心建立在对模型每一个细节的有力把握。
  14. 量化投资一个人是搞不定的,哪怕你是明星基金经理。
  15. 要有异常值的处理机制,不要小看异常值。
  16. 不要轻易相信别人做出来的结果,条件许可的话最好自己亲自动手做一做。(2012-02-09)
  17. 一个小错误可能导致致命的结果,就因为多添加了一列数据。(2012-02-29)
  18. 模型用不同的数据结构也许产生相同的结果,但是效率却可能大相径庭。比如考虑两个因素的时候,用矩阵可能比较好;而考虑三个及以上因素的时候,数据框可能更快捷。(2012-05-07)

因死而完美

昨夜看《平凡的世界》有点伤心,田晓霞突然被死掉了,一时让我难以接受,而孙少平比我更难以面对这个事实。田福军两眼含满泪水对孙少平说:“孩子,……,在她活着的时候,你曾给过她爱情的满足……我深深地感激你!她留给我们的主要纪念就是十几本日记。这三本是记述你们之间感情的,就由你去保存……”

也许,精神上两相执手而世俗中难以终成眷属的恋人,一方死去未尝不是一种永恒相爱的完美结局。

看到这里,感觉《平凡的世界》是个悲剧,但是比喜剧更荡气回肠、令人动容;而比之“伟人”,孙少平这样的“平凡人”更让我起敬与共鸣。

中国统计之网

有些时间没写学术的东西了,老是儿女情长、无病呻吟,感觉还是不靠谱。

现在的工作其实有点学术,但是涉及商业机密,所以不便写出来,暂时发布一点曾经研究过的玩意儿。这篇博文的内容其实在之前的《统计人话》中蜻蜓点水地提到过,由于当时处于毕业季,怕惹麻烦,所以一直没有真正公开,不过在统计之都内部邮件中交流过,老谢说“这篇文章对COS肯定有很大的宣传作用”,但我一直拖着没发,实在对不住老大。现在不妨在这里先试一下水吧。言归正传,开始学术:

一、数据来源与处理

下面将要分析的数据是从中国知网抓取的中文统计核心期刊《中国统计》(2000.1 – 2010.11)、《统计研究》(2003.1 – 2011.1)以及《数理统计与管理》(2000.2 – 2011.2)十年左右的数据。原始数据都存储在Excel表格中,其格式如下:

表1中从第二行开始,每一行是一篇文章的记录,上述三大期刊分别有729、1897、1371条记录,一共3997条记录。

这里要研究的是合著关系网络,因此需要将表1中的第二列的作者数据提取出来。不难发现,表1中当一篇文章的作者多于一个的时候,是使用分号“;”隔开的。我们可以借助这个分隔符进行中分分词,这样就能将所有作者提取出来。在搜索每篇文章以提取作者的时候,顺便去掉了没有作者的文章(这种文章一般不是论文),去掉了只有一个作者名但其字数大于4的文章(一般是以某个机构或小组署名的文章),以及去掉了英文名。于是得到了如下形式的每篇论文的作者列表(一共3617篇文章,以前6篇文章为例):

经过拼音排序后没有去重的作者列表(一共6568位,以前20位为例):

不同作者发表论文的频次表(一共4293位不同的作者,以后50位为例):

二、基本统计与网络分析

 2.1 描述统计分析

表5是上述中文统计核心期刊的一些基本统计数据:

由表5可以看出,平均每位作者发表1~2篇论文,极个别人发表了33篇论文。而从图1中我们可以发现,在核心期刊发表1~3篇文章的人占了绝大多数,能发5篇以上的人就非常少了。

由表5和图2可以看出,一个人发文章或者两到三个人合写文章的情况比较常见,一篇文章的作者数超过4人的情况极少。值的注意的是,有一篇文章的作者数达到了18人,通过查找原始数据,发现它是一篇国际统计论坛的综述,其作者是该论坛主办单位的许多老师。

我们还能从表5和图3中发现,平均每位作者大约有1~2位合著者,“单干”的情况也较常见,为数不多的人会有10个以上的合著者。

通过这些图、表的简单展示和分析,我们对中国统计学界近10年的情况有了一些粗浅的了解。

 2.2 权的定义与合著关系网络的构建

下面我们来讨论中文统计核心期刊合著关系网络的构建问题。

合著关系网中的顶点是全部论文的所有作者。如果两位作者的名字同时出现在某一篇论文的署名当中,那么就在两者之间连接一条边,表示他(她)俩有过合著关系。如果一篇论文有n个作者,那么就产生了n(n – 1) / 2种两两合著关系。比如表1中的第2篇论文就可以建立3(3 – 1 ) / 2 = 3条边:周建章–余洁平,周建章–苏毓腾,余洁平–苏毓腾。注意,当两个特定的作者vw同时出现在多篇文章当中时(这是俩人合著关系强弱的某种体现),我们对vw不再重复连接,这样保证了网络是简单网络。

如果我们只关心作者两两间有没有合著过,而不关心合著关系的强弱,那么我们得到的网络就是一个无权网络。但是权重信息对于发现有效的社群结构可能是非常有用的。所以我们有必要将合著关系的强弱体现在边的权重上,从而构造一个加权图来更好地实施合著关系网的社群挖掘。

接下来的问题就是:怎么构造边的权重?

曾经有人利用两个作者合著的次数来定义权重。这种定义方式考虑了合著次数这个反映关系强弱的重要信息,但是它有不足之处:根据我们的常识,一篇两个作者的论文很可能比一篇6个作者的论文有更强的两两合著关系,因为作者之间的熟悉程度与时间投入一般会因为作者数目的增多而减弱。

我们综合考虑合著次数和单篇论文作者数目这两个因素来定义作者vw之间权重:

其中nk表示论文k的作者数目,当作者v在论文k中出现时,\delta_v^k等于1;否则为0。 公式中的分母nk – 1排除了只有1个作者的情况(即v不能等于w),因为这种情况对合著关系来说是没有意义的。

图4是该定义的一个示例,图中的作者vw一起发表过3篇论文,其中第一篇论文有4个作者,第二篇论文有2个作者,第三篇论文有3个作者,于是vw在3篇论文中的合著关系强度分别是1 / 3、1、1 / 2,所以vw总的合著关系强度是1 / 3 + 1 + 1 / 2 = 11 / 6。

我们计算出最大权重是8,最小权重是0.06,可见合著关系最强与最弱相差悬殊。而平均权重是0.65。

定义了顶点、边与权重,那么加权合著关系网络的构建也就完成了。我们构建的网络有4293个顶点,3827条边。而理论上最多有(42932 – 4293) / 2 = 9212778条边,可见这种合著关系网络是很稀疏的,其实稀疏性在社会关系网络中是很常见的。

2.3 聚类系数分析

聚类系数分为局部聚类系数和全局聚类系数,它们反映了网络中局部或者全局的群聚效应。

我们计算出合著关系网络的全局聚类系数是0.57,跟Newman文献中的7种合著关系网的聚类系数相比处于第2位(第1位是0.726,第3位是0.496)。可见从整体上来说,中国统计学界合著的群聚效应比较明显。

我们还发现,网络中有2396位(约占整个网络的56%)作者的局部聚类系数不存在,即他们的度为1,也就是说他们都只与1个人合写过文章。令人惊讶的是,有1359位(约占整个网络的32%)作者的局部聚类系数为最大值1,可见网络中存在相当多的“三角”合著关系。有151个人的局部聚类系数为最小值0,也就是说与这些人合著过的作者之间并没有发生合著。

2.4 距离与中心性分析

2.4.1 距离与最短路径

两个顶点间的距离就是它们之间的最短路径的长度。现在我们考虑的是加权图,所以距离不再是边的条数之和了,而是边的权重的倒数之和。之所以取倒数,是因为权越大其倒数就越小,意味着关系就越近。

我们计算出最大加权距离为42,平均加权距离为16。这种距离虽然考虑了权重,但是不够直观。反而无权距离的实际意义更明显一些,无权距离减去1就表示两个人之间建立联系至少需要经过多少个人。我们算出来的无权距离最大值是18,平均值是6.8,后者表示在该合著网络中任意两位作者建立联系大约平均至少需要经过5~6位作者,这与社会学中有名的“六度分离”理论比较吻合,表现出“小世界效应”,这也许反映了合著关系网络一定的社会性。

2.4.2  中心性

中心性是社会网络中的一个重要概念,它是度量网络中成员重要性的某种指标。常见的中心性包括:顶点的度、中间性、接近性等等。这些中心性指标从不同的角度反映了顶点的某种重要性。

顶点的度在无权网络中是与该顶点连接的边的条数。在加权网络中对应的概念是顶点的强度,即与该顶点连接的边的权重之和。有意思的是,在合著关系网中,顶点的强度其实是作者合写文章的数目:

我们的合著关系网中,顶点强度最大值是20,最小值是0,平均数是1.2,中位数是1。顶点的强度排名前10的作者如下表所示:

顶点的中间性是经过该顶点的最短路径的条数(给定顶点对之间若有多条最短路径,则不重复计算)。顶点的中间性反映了顶点起“桥梁”作用的大小。不过我们现在将最短路径替换为加权最短路径。顶点的中间性排名前10的作者如下表所示:

顶点的接近性是某顶点到所有其他顶点的平均最短路径长度的倒数:

其中| V |是网络中总的顶点数,gd(v, w)表示顶点vw的最短路径的长度。cnv的值越大表示顶点v与网络中的其他成员越接近。我们的加权合著关系网络中顶点的接近性的最大值为110.9,最小值为0(这是孤立点),平均数为7.4,中位数为1.5。平均数比中位数大很多,意味着接近性呈偏态分布,大多数人的接近性较低,少数人的接近性很高。顶点的接近性排名前10的作者如下表所示:

比较表6、表7和表8,我们会发现,在3张表中都出现过的人只有谢邦昌(台湾辅仁大学统计信息学系教授)、朱建平(厦门大学计划统计系主任)和曾五一(厦门大学经济学院副院长、中国统计学会副会长);出现过两次的只有缪柏其(中国科学技术大学前任统计与金融系主任)、金勇进(中国人民大学统计学院前任院长)。数据表明,这些学者不但中心性高,而且发文数也高,而实际上他们也确实是中国统计界很有影响力的人物。尤其引人注目的是,谢邦昌教授的3种中心性都排名第一,即他与别人合写的文章最多、“桥梁”作用最大、与同行最“接近”。需要说明的是,由于我们只选取了3个中文核心期刊10年左右的数据,还没有考虑英文期刊,所以很可能会遗漏一些重要的学者以及与真实情形会有偏倚。

值得注意的是,表6中顶点的度的高低与发文数是正相关的,但其他表中的中心性就不全是这样的了。有些中心性高的作者发文数并不多,比如李海涛(曾五一的博士研究生)和来升强(朱建平的博士研究生),还有李善同(国务院发展研究中心发展战略和区域经济研究部研究员)。李海涛和来升强的导师都是中心性很高的人,所以可能对二者的中心性有很大的拉动作用;而李善同与表7中的袁卫(中国人民大学前任统计系主任)和表6中的许宪春(国家统计局副局长)都合发过文章,她的中间性很可能因此提升,而且她对统计学与经济学的交叉研究起到了较好的桥梁作用。

综上我们可以归纳一点:多发文章或者多与别人(特别是有影响力的人)合写文章,是在学术界中建立广泛联系并提升自身影响力的一种好方式。

 三、社群挖掘分析

3.1 合著关系网中的成分

前面我们对加权合著关系网络进行了探索性描述分析和传统社会网络分析,获得了对这个网络的一些基本的感觉与认识。下面我们尝试利用WCMN算法(加权CMN算法)对该合著关系网络的社群结构进行挖掘分析,以期得到有意义、有意思的学术社群结构。

我们对社群有一个基本假设,就是社群内部的每一对顶点之间至少存在一条内含的路径,也就是社群要有连通性。我们主要对非连通图的每一个连通成分进行单独分析,从而简化了在非连通图上挖掘社群的复杂性。一个连通成分就是一个局部最大的连通子图。

我们的合著关系网络有4293个顶点,3872条边。经检验,该网络不是连通的网络。这是很正常的,这么大的社会网络中任意两个成员之间都连通几乎是不可能的。但是我们的网络中存在1773个连通成分,具体情况下表所示:

从表9中我们可以看出,绝大多数连通成分包含的顶点数都只有区区的几个。唯有两个连通成分的规模明显大于其他连通成分:一个含有383个顶点,另一个含有156个顶点。

我们的第一大连通成分包含383个顶点,750条边,如图5所示。下文我们将以最大连通成分为例进行社群挖掘。在进行社群挖掘之前,我们从这些网络中几乎是看不出什么社群结构的。

3.2 模块性与最佳社群数

WCMN算法是一种直接对模块性Q值进行最优化以寻找最佳社群结构的凝聚算法,它朝着使Q值增大最快或减少最慢的方向将社群一步一步地融合。在这个过程当中会出现Q的峰值Qmax,如图7所示。我们认为Qmax对应的社群结构就是最佳社群结构, Qmax越大表明我们的算法效果越好。


上图是WCMN算法对最大连通成分进行社群挖掘的Q值变化图。因为有383个顶点,所以一共有382次融合。在第359次(即倒数第24次)融合的时候Q达到最大值0.87,所以最佳社群数目是24。这里的Qmax高达0.87,所以这个最大连通成分中可能存在明显的社群结构,同时也可能是WCMN算法效果良好的一个表现。

3.3 合著社群及其可视化

我们已经利用WCMN算法将最大连通成分的383位作者划分为24个合著社群,如表10所示,表中每位作者上方的数字是1~383中的某个编号,编号大小与作者在我们所选的期刊中的发文数成正相关。通过2.4节的分析,我们知道这里编号较大的作者一般处于网络的中心位置或者桥梁位置,观察下文的图8发现确实如此。

为了更直观地了解这24个社群,我们将383位作者构成的网络进行可视化展示,如图8所示。图中一共有24种渐变的颜色,每种颜色代表一个合著社群,圆圈的大小与顶点的强度正相关,线条的粗细与边的权重正相关。

图8   WCMN算法社群挖掘效果可视化

分析表10和图8,可以发现我们挖掘出来的社群与实际情形比较吻合,这些社群里面蕴藏了很多有意思、有价值的信息:

  1. 反映了不同的学校或研究机构:中国人民大学统计学院(社群1),国家统计局(社群2),中国科学院数学与系统科学研究院(社群3),吉林大学(社群6),山西财经大学统计学院(社群11),厦门大学经济学院计划统计系(社群14),中国科学技术大学统计与金融系(社群15),西南财经大学统计学院(社群16),东北财经大学统计系(社群17),暨南大学经济学院统计系(社群19、20),北京航空航天大学经济管理学院(社群21),北京大学光华管理学院(社群22),清华大学房地产研究所(社群24)等。
  2. 反映了不同的研究方向或领域:国民经济核算与投入产出分析(社群2、23),国际竞争力(社群8),宏观经济统计(社群10、11、14、17),金融统计(社群3、15、21),数量经济(社群6),质量管理(社群4、7),数据挖掘(社群9),风险管理与精算(社群12),数理统计(社群13),抽样调查(社群18),社会经济调查(社群19),绿色核算(社群22),房地产经济(社群24)等。
  3. 反映了不同机构、学者或领域之间的沟通与合作:国家统计局与高校的合作(社群2、4),台湾、人大、厦大、清华在数据挖掘领域的合作(社群9),人大与天津财大在风向管理与精算方面的合作(社群12),不同单位在质量管理方面的合作(社群7),不同机构的数理统计学者的合作(社群13),数学学院与商学院的合作(社群6),数学与金融的结合(社群3),抽样技术与数据挖掘的结合(社群9),房地产与统计的结合(社群24),老师与学生的合作(社群8、15、18、19、20、23),学者与母校的合作(社群13、16)等。

以上结果很好地印证了前面的说法:高达0.87的Qmax值预示着这个网络中可能存在明显的社群结构。同时也有力表明:WCMN算法对这种加权合著关系网络具有良好的应用效果。

如果不考虑权重的话,那么WCMN算法就变成了CMN算法,这时候我们得到Qmax = 0.79,最佳社群数是20,不出所料两个数值都比前者的小。效果其实也还不错(如图9所示),只是没有WCMN算法挖掘得那么精细,毕竟它利用的信息少一些,比如人大统计学院的国际竞争力、抽样技术、风向管理与精算三个研究方向就没有区分开;还有将谢邦昌老师归于人大统计学院这一群体似有不妥,而WCMN算法是将谢老师与人大、厦大的某些数据挖掘方向的老师与学生整合为一类。

图9   CMN算法社群挖掘效果可视化