这是一个豆瓣上的小站，现在还处于早期阶段，我们希望把它建设成为一个国内对互联网算法感兴趣的同学们（简称算法攻城师）交流的基地。现在暂时由文栋、xlvector、yoyo和我打理，欢迎大家加入。

由于处于建设初期，版块刚刚确立起来，内容还不甚完善，将来会有更多的内容更新上去，也欢迎大家在上面交流心得，讨论问题。

协同过滤版块，会转载一些范畴广泛的涉及互联网算法的博客文章，可由大家推荐产生，转载过来（保留回链），然后大家就该文进行交流，这样也方便把不同地方的精华内容汇聚到一起进行沟通，省了大家东奔西走之苦。

新产品新应用、学术动态、电子杂志、招聘广场各版块，顾名思义，都有其各自的职能，不过现在还有待开发，希望将来有会不同的版主来管理，大家可以积极申请。

小站的首页会用于将来线下活动的召集，活动反馈，讨论，内容建议等等。

希望与大家一同把这个基地建设好。

举个简单的例子来说明，关联规则一个经典的应用是计算超市购物中被共同购买的商品，它把每个顾客的一次交易视作一个transaction，计算在不同transaction中不同item组合的规律性。而如果我们考虑一个用户多次在超市购物的情况，那么这些不同时间点的交易记录就构成了一个购买序列，N个用户的购买序列就组成一个规模为N的序列数据集。考虑这些时间上的因素之后，我们就能得到一些比关联规则更有价值的规律，比如关联挖掘经常能挖掘出如啤酒和尿布的搭配规律，而序列模式挖掘则能挖掘出诸如《育儿指南》->婴儿车这样带有一定因果性质的规律。所以，序列模式挖掘比关联挖掘能得到更深刻的知识。

在实际当中，序列模式挖掘被广泛地应用于各种序列数据集中，如生物信息学上的基因微阵列数据，从中挖掘哪些基因组合模式在某类病人中会频繁出现；以单词作为item的文档序列，研究在不同文档中单词序列的出现模式；用户点击流数据，用于挖掘用户的频繁点击模式，建立用户模型，完善网站功能与UI结构。除此之外还有很多，只要是序列数据集，都可以考虑利用序列模式挖掘获得规律。

图1是一个序列数据库，及其以0.75作为最小阈值（min_sup）的频繁序列模式。借此介绍序列挖掘中的几个主要概念。

图1

• 序列（Sequence）：以SID表示，一个序列即是一个完整的信息流。
• 项目（Item）：序列中最小组成单位的集合，比如在这个样例中的项目为{A, B, C}。
• 事件（Event）：通常用时间戳标志，标识事件之间的前后关系。又叫Itemset，是Item的集合，样例中以EID表示。
• k频繁序列：如果频繁序列的项目个数为k，则称之为k频繁序列，以Fk表示（图1的F1,F2,F3）。
• 序列的包含关系：对于序列x和y，如果存在着一个保序的映射，使得x中的每个事件都被包含于y中的某个事件，则称为x被包含于y（x是y的子序列），例如序列B->AC是序列AB->E->ACD的子序列。
• 支持度（support）：某序列x的支持度是指在整个序列集中包含x的序列的频次。

有了以上概念之后，序列模式挖掘的问题就定义为：给定一个序列数据库以及最小支持度min_sup，找出所有支持度大于min_sup的序列模式。

解决这个问题的一个著名的算法叫SPADE，以及加入各种约束限制后的cSPADE，它们来自同一个作者Zaki。SPADE是无约束的序列频繁模式挖掘算法，其基本思想是利用一个S后缀等价类（suffix equivalence class）的概念[S]，从F1开始计算，由F1张出以F1中的频繁序列为[S]的F2，再以F2张出以F2中的频繁序列为[S]的F3，以此类推。由于引入了后缀等价类这个概念，使得整个算法对数据库的遍历次数与计算量大大降低，是一个现实中有效的算法。cSPACE在前者的基础上引入了对序列模式的约束，其中包括：

• 序列长度与宽度的约束（序列的长度是指序列中事件的个数，宽度是指最长事件的长度）；
• 最小间隔的约束，即事件之间的最小时间间隔；
• 最大间隔的约束，即事件之间的最大时间间隔；
• 时间窗的约束，即整个序列都必须发生在某个时间窗内；
• 其它，如只限制为某些Item，序列数据库存在分类标签。

以上的约束基本已经覆盖了实际应用中所有可能碰到的限制情况，实际上约束并不常用到，通常SPADE算法已经够用了。

由于数学表述过多，而且整个过程比较繁琐，这里不打算详细叙述两个算法的计算流程，有兴趣的可以参看文后附上的参考文献。

很幸运的是，万能的R已经有序列挖掘的包，对SPADE算法与cSPADE算法有了完整的实现，即使对算法并不充分了解，也可以把它们用于自己的场景中。它就是arulesSequences，以下以一个例子简单介绍一下它的用法。

1
2
3
4
5
6
7

library(arulesSequences)
 data(zaki)
 s1 <- cspade(zaki, parameter = list(support = 0.4), control   = list(verbose = TRUE))
 summary(s1)
 as(s1, "data.frame")
 s2 <- cspade(zaki, parameter = list(support = 0.4, maxwin = 5))
 as(s2, "data.frame")

第2行导入数据库zaki，它的结构类似于图1所示的那种序列-事件集，第3行设定最小支持度为0.4，并用无约束的SPADE算法搜索频繁序列模式，第6行引入了时间窗的限制，用cSPADE算法进行搜索。整个过程十分简洁明了。

更多使用的例子可在R的终端用?cspade命令查看。

参考文献：
Sequence Mining in Categorical Domains: Incorporating Constraints, Mohammed J. Zaki, 2000 ACM

一个人的旅程，总会有诸多的不便，却也总能发现一些额外的风景，比如旅途中的人，旅途中的感觉。因为孤独，所以更容易与别人沟通来。因为安静，更容易有所思。

车上的老乡

北京到丰宁的漫长车程之后，丰宁到坝上的路程仍然很难熬，于是车上我认识了一位当地的老乡，虽然仍旧神采奕奕，但其实已经七十出头了。他不断地跟我说着他年轻时骑马一天奔跑几十公里去另一个地方的故事，又给我传授骑马的要领，以及坝上这几十年来的变迁。岁月的沧桑留给他额头上深深的皱纹，但从他豪迈的言谈当中，我仍能依稀地找到那个当年在马背上意气风发的少年。有一天我也会老去，回想起年轻的这段岁月，我是否也会像他那般充满豪情，充满骄傲。这也许就是我正在追求的。

那一家子

他们是我到坝上后结交的第一批游客，我们住在同一个农家院里。当我下午四点拖着疲惫的身体和心理来到住处时，他们刚刚吃完饭，一对夫妇，一个七八岁的孩子，一个六旬老人。当我听说他们在之前一个晚上已经出发，却因为路上塞车而只比我早到了半个小时时，我心里稍稍找到了些许平衡。天公不作美，长途跋涉的我们遭遇的不是一个明亮的草原，而是一个雨水纷飞而不得外出的阴霾，共同的行程共同的遭遇让我们很快便聊到了一块，从来程计划，到归程计划，没有新意却饶有趣味。第二天的一大早，为了弥补前天的遗憾，我起了个大早，却发现他们起得比我更早，因为我们都需要在中午吃过饭后便去赶车，所以在争分夺秒这一点上我们又是一致的。早餐前趁着晨光和他们一道在村边山上游玩，本来还约定早饭后要一起租向导一起骑马一起坐车回北京的，不料我座上那任性的青马一出村而脱缰狂奔，把我带到了另一个队伍，短暂的缘份只能就此终结了。我的相机里，还保存着他们一家四口的照片，兴奋写在了每一个人的脸上，实际上，他们的情绪并不单体现在镜头前，他们的情绪也感染了我。我一直都没想明白，到底是什么让如此奔波的他们毫无怨言反而精神百倍地生活着，是广阔草原的魅力，是出游的兴奋，还是家庭的温馨相伴？我体验不到，因为我只是一个人而已。

两兄弟

周六的晚上，农家院为住客准备了烤羊腿。一个中年男子熟练地拿刀子和酱料给炭火上的羊腿做着细致的功夫。第一次面对这种情景的我颇有些手足无措，他就很热情地用刀子帮我割出表面已经烤熟的肉，分给我吃。如果不是后来与他的弟弟聊起天，我肯定会一直认为他就是农家院的主人。事实很让我吃惊，由于他每年都会来这里几次，所以他已经把这里当成家一般了，每次都操刀烧烤，宛如家中的主人。弟弟也已近中年，极其健谈，坐在我空荡荡的房间里就跟我侃起了大山，从海南说到山海关，又从上海说到内蒙古。阅历甚广，而且事业有成的他经常会带上一家大小，享受着旅行的乐趣。这次本来并不太想来，但为了这个哥哥，也就陪着走了一趟，即便不出去玩，看着家中老小各得其乐，也便就可以了。哥哥处处可为家，弟弟走到处处心中都有家，这也是我很难体会到的，因为离家十七年来，我再也未曾有过家的感觉。漂泊有定数，归家无可期。

同骑共行

实际上最后我并没有去赶回北京的长途车，因为遇上了两位自驾的游客，并坐了他们的顺风车回京。本来我不太可能会跟他们结识的，虽然我们也住在同一个农家院里，但是我那脱缰的青马却准确地把我带进了他们的马队里，于是便有了后面的故事，这就是所谓的缘份吧。当三匹快马狂奔在草原上，穿过跑马场，穿过鲜花烂漫的草地，跑到白桦林，停留在中途驿站喝茶，我们自然而然地多了很多的共同话题。也凭着这份萍水之交，他们邀请我坐上了他们回京的车。我后来想，凭个性而论，如果在城里，我也许不太可能与他们结交，但在一片远离大都市的土地里，在一个陌生的地方，人们似乎更容易建筑起感情来。在广阔的天地之间，我们没有了芥蒂，所有的人都只是天地渺小的儿女，我们都需要伙伴，需要温暖。

八年前长江边上，我与那几个共船共渡的朋友便是如此。我至今仍清晰地记得那山西大哥及他的妻女，夕阳下长江边他们的嬉戏，以及宜昌火车站旁住满了人的小旅馆，半夜里，他起床帮我们买票，而我守着他的家人。那个叫王硕的小女孩，现今应该也长大成人了。

大部分的人注定都只是你生命中的匆匆过客，但有一些会成为你记忆中永恒的风景。

上述对信息的解决方案都是通过follow人来实现，而follow内容无疑是另一条可行的路径。关于内容的分类组织已经有很多年的研究与实践，在郑昀的这篇关于Topic Engine的博客中有很详细的综述，而对这些分类内容的follow，以得到一个类别的信息更新，就我所知，还并不多见。比较知名的如google资讯与google快讯，前者属于对内容的分类组织，后者则属于对分类内容的订阅或者说follow。依托于google强大的搜索能力，信息的新颖性及广阔性很有保证，但由于产品的定位并非要作一个详细的内容分类，所以分类比较粗糙，只是一些如门户网站般的粗分类别而已。

没有了四年前的激情澎湃，这一次，我大部分的场次，都是坐在那个10cm*15cm的网络直播窗口前，听着那一成不变的苍蝇声，看着时间静静地在球场每个角落流淌，看着岁月的痕迹划过卡纳瓦罗的眼眶亨利的额头，看着成熟与责任写满了C罗与梅西曾经青涩的脸庞。我们也在随着他们一个四年又一个四年地长大，变老。

我必须得很俗地冒用那句流行的网络用语：世界杯就是个大茶几，个别人得到洗具，大多数人只能得到杯具。让我们简单地回顾一下，这一个月的时间里，他们所演绎过的人生。这些人，存在于世界的每个角落，这些事，正在我们的身边发生。

亨利：即使奉献所有，组织也会毫不犹豫地把你抛弃

如果这个世界上有如果，亨利还会不会选择冒着搭上自己一世英名的风险，把法国队送到南非？他是个天才，他也很勤奋，只是他一直都生活在巨人的阴影下。直到某一天他认为自己可以当家作主了，才忽然发现，自己已经老了，老得即使你为集体牺牲了所有，集体也会毫不犹豫地把你牺牲掉。而自己忍辱负重所换来的，只是另一个更耻辱的结果，最可悲的是你只能眼睁睁地看着它发生，除了做一些场下的努力，然后被所有的人唾骂，你无能为力。前些天，亨利宣布退出国家队，同多梅内克一样，他的错误只是在于没有在合适的时间选择离开，同样情况的还有卡纳瓦罗，我们不能责怪他们，就在大家都在赞颂急流勇退的勇气的时候，他们心中只是仍有一些信念在坚持，这在足球场上很常见，一如大卫贝克汉姆。

科威尔：徘徊在宿命与挚爱的人

本来一次无可争议的判罚，却因为这个人背后的故事，而变得无比悲情，他可以坦然接受命运的判决，却无法释然地接受裁判的判罚。这是本届世界杯最触动我的一个人、一件事，我看到了一个科威尔身后千千万万个科威尔的身影。疾病几乎夺走了他的一切，除了他的追求，所以他的理想得以继续，直到他再一次倒在命运的嘲弄下，所谓造物弄人。为了保护我们脆弱的心灵，故事、电影，常会给我们美好的结局，然而这才是现实，他抗争了，但是他再一次失败了，如果他选择继续抗争，他也许还会继续失败下去，直到有一天，他松开那根死死缠着他的绳索。

马拉多纳：感性与理性的交集在哪里？

哪一名球员是世界杯上最耀眼的明星，这个问题可以众说纷纭，但如果问世界杯上谁最具有感染力，那毫无疑问就是马拉多纳。这个一生充满争议的天使与魔鬼的结合体，在一次0：4的惨败之后，也再一次深深地跌入到自己情绪的谷底。疯狂如马拉多纳，即使他失败，他的队员仍然愿意为他而继续战斗，他的国民仍然选择相信他。只是这一回，是彻底地倒在了理性、科学的德国人脚下，再无如二十年前那般争议的余地。当性感的探戈遭遇严肃的战车，失败的总是感性的舞者，更何况，现在这个舞者根本都无法轻盈地舞动。即使阿根廷人仍然会选择相信这个他们心中神一样的男人，相信他必定会找到一个感性与理性的交汇点，然而，二十多年前，阿根廷依靠的是马拉多纳，二十多年后，阿根廷依靠的还是马拉多纳，阿根廷足球和巴西足球的差距，也许就是一个不世出的天才与很多很多天才的差别。

苏亚雷斯：当童话发生在人间

比赛的最后一分钟，作为一名前锋的他，作出了排队运动员般的拦网动作，把必进球扑了出去，又在沮丧中目睹了自己球队从地狱到天堂的过程，自己则经历了从罪人到英雄的一幕，经此戏剧性的几分钟，谁能不动容。世界上无数人被这一幕感动了（非洲人民应该除外），这一刻，人们忘记了身边一切应当有的秩序，而只有时间去感叹：这是一件多么奇妙的事情啊！我想，这其实就是童话故事在有条不紊规规矩矩的世界中的呈现，就像多年前的黛安娜王妃，当灰姑娘出现在自己的世界，谁又能不感动。当丹麦童话、希腊神话在欧洲杯赛场上上演，你会知道，在足球场上，童话是会真真切切地发生的。那场比赛结束后，我曾发了一条广播调侃道：按照近年南美和欧洲轮流夺冠的规律，我猜乌拉圭进决赛，苏亚雷斯在决赛中以守门员身份出战，最后点球大战扑出对方5个点球，然后自己罚进最后一个球，成为一代传奇。奇迹真的会有的，童话真的会发生的，如果在你的生活中没有，那么，关注那个绿茵场吧。

巴拉克：你以为自己是属于她的，她却渐渐远离了你

虽然他没有来，虽然在足总杯的决赛中他被自己德国队的队友博阿腾的哥哥加纳的博阿腾铲伤而无法参加这也许是他最后的一次世界杯，但在很多人眼中，这些年的德国队与他已经难以分清，从他伤后到德国训练营探营，从他到南非看德国队的比赛，从拉姆说即使他回来自己也不会交出队长袖标，从他在半决赛前因为倍感冷落离开南非并坚持声称如果德国进了决赛还会回来，从德国结束所有比赛后施魏因施泰格说自己心中的队长只有巴拉克，谁又能区分德国队和巴拉克。随着德国足球风格从传统坚韧硬朗到现在华丽奔放的转变，他很有可能成为旧德国战车的最后一个代表性人物。当德阿大战中德国新13号托马斯穆勒大放异彩，解说员说如果镜头没有捕捉到巴拉克坐在观众席上，大家也许早已经忘记这位曾经的13号的主人。从来只见新人笑，哪里见得旧人哭，切尔西已经不再需要他，德国队也已经不再需要他，他一直在追求自己的巅峰时刻，却永远都只差一步，最后剩下来陪伴他的，也只是落寞。真的斗士，即使落寞如是，也会选择坚持下去。

又一个四年落下帷幕，又一个四年将揭开它神秘的面纱。面对着过去的和即将来到的一切，我们又将何去何从。我想，三四名比赛第92分钟所发生的就浓缩了世界杯关于结束的一切，当弗兰的任意球呼啸着击中横梁，裁判的哨声随即响起，在高潮中落幕，在紧张中分出了胜负，留给你无尽的想像，而你还未曾有时间从中走出来。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

# 读入数据，原数据是user-subject的收藏二元组
data = read.table('data.dat', sep=',', header=TRUE)
# 标识user与subject的索引
user = unique(data$user_id)
subject = unique(data$subject_id)
uidx = match(data$user_id, user)
iidx = match(data$subject_id, subject)
# 从二元组构造收藏矩阵
M = matrix(0, length(user), length(subject))
i = cbind(uidx, iidx)
M[i] = 1
# 对列向量（subject向量）进行标准化，%*%为矩阵乘法
mod = colSums(M^2)^0.5  # 各列的模
MM = M %*% diag(1/mod)  # M乘以由1/mod组成的对角阵，实质是各列除以该列的模
#crossprod实现MM的转置乘以MM，这里用于计算列向量的内积，S为subject的相似度矩阵
S = crossprod(MM)
# user-subject推荐的分值
R = M %*% S
R = apply(R, 1, FUN=sort, decreasing=TRUE, index.return=TRUE)
k = 5
# 取出前5个分值最大的subject
res = lapply(R, FUN=function(r)return(subject[r$ix[1:k]]))
# 输出数据
write.table(paste(user, res, sep=':'), file='result.dat', quote=FALSE, row.name=FALSE, col.name=FALSE)

除去注释，有效代码只有16行。其中大量运用了向量化的函数与处理方式，所以没有任何的显式循环结构，关于向量化更详细的叙述可看这里。

注：该代码实现的只是最基本算法，仅作参考，不承诺在大规模与复杂数据环境下的实用性。

源数据文件data.dat的内容如下所列：

user_id,subject_id

1,1

1,3

1,7

1,13

2,2

2,5

2,6

2,7

2,9

2,10

2,11

3,1

3,2

3,3

3,4

3,7

3,9

3,10

5,13

6,1

6,3

6,4

6,5

6,8

6,10

8,1

8,2

8,3

8,5

8,6

8,7

8,8

9,13

10,12

11,2

11,3

11,4

11,6

11,8

11,9

11,13

12,12

13,3

13,6

13,7

15,4

15,12

15,13

16,2

16,3

16,4

16,7

16,8

17,2

17,3

17,4

17,5

17,6

17,7

17,8

17,9

17,10

17,11

18,2

18,3

19,2

19,3

19,5

19,6

19,9

19,10

19,11

19,12

20,1

20,3

20,4

20,7

20,13

21,1

21,6

21,8

21,9

21,11

21,12

21,13

22,6

23,2

23,4

23,9

23,12

24,1

24,5

24,9

25,2

25,6

25,10

25,11

26,2

26,3

26,8

27,3

27,6

27,12

27,13

28,1

28,2

28,3

28,5

28,7

28,9

28,10

28,11

28,12

28,13

29,1

29,2

29,3

29,4

29,5

29,6

29,7

29,8

29,9

29,10

30,6

30,7

30,9

30,13

31,6

31,11

32,1

32,5

33,2

33,13

34,3

34,7

34,8

34,9

34,10

34,13

35,3

35,4

35,5

35,6

35,7

36,2

36,3

36,4

36,6

36,7

36,8

36,9

36,11

36,12

36,13

38,5

41,1

41,3

41,4

41,5

41,6

41,7

41,11

42,2

42,3

42,7

42,8

42,9

42,10

42,11

43,2

43,6

43,10

43,11

43,12

我最初接触到专门用于科学计算的工具，是大名鼎鼎的matlab，正如它帮助了无数中国学生顺利毕业的赫赫功劳一样，它是我对于向量化计算的启蒙老师。用过matlab的人都会对其循环结构的效率无法忍受，不知道是否有意而为之的这样的设计缺陷，迫使人们要想真正地用好它，就得接受它提供的向量化计算的思想，在掌握了这个专门为高效计算而设计的计算思想之后，你会发现自己获得的不单是计算效率上的极大提升，还有算法设计思想上高屋建瓴式的跨越。

毕竟matlab是更适合于研究领域的商业软件，而且是闭源的，毕业后不久，我就选择了R作为matlab的替代品，看中的正是R在向量化计算支持之余的灵活性及丰富的第三方库，似乎天生就是数据挖掘人员最趁手的那把刀。因为我需要的就是这样的一个计算环境，它不单是一门编程语言，也不必是一个已经很完备的工具。它就是这样的一个环境，拥有很多的各领域相关的工具包，我可以不必操心太多过于底层的或与工作主题没有直接关系的问题，同时当不满意时我又具有对它最自由的掌控。实际上，我寻求的就是matlab与C之间的一个平衡点。R是一个面向科学工程计算特别是统计计算的工具，与matlab一样，其循环结构的效率也无法让人满意，所以，我们必须利用向量化的编程范式，必要时采用并行/分布计算（因为，向量化本质上就是一种并行计算，也是我们通常理解那种并行计算的天然先驱）。而这一切，在R面前，都不是什么问题。

所谓向量化：

wikipedia中的定义是：Vectorization is the more limited process of converting a computer program from a scalar implementation, which processes a single pair of operands at a time, to a vector implementation which processes one operation on multiple pairs of operands at once。向量化计算是一种特殊的并行计算的方式，相比于一般程序在同一时间只执行一个操作的方式，它可以在同一时间执行多次操作，通常是对不同的数据执行同样的一个或一批指令，或者说把指令应用于一个数组/向量。

像R和matlab这样的现代科学软件包，都会以向量化作为其计算的基本特点（即使python的numpy包也是如此），所以在R的基本运算中，随处可见向量化计算的影子，以下仅举几例，以让读者了解向量化是多么地简单和直观：

向量取值：V[1:10]  （把get操作应用于向量V的不同元素）

向量赋值：V[1:10] <- seq(1,10)  （把set操作应用于一个序列与向量V的对应元素）

apply系列：lapply(V, mean)  （跟python的map函数类似，是向量化最直接的表达形式）

矩阵运算：A + B；A %*% B  （矩阵的基本运算也是向量化的典型形式）

这些操作很常见，以致于我们都没有意识到这就是一种有助于提高计算性能的向量化计算，更忘了由此而把这样的思想扩展到我们算法设计的更多方面。熟练使用它，你获得不单是计算上的好处，还有对问题理解上的整体性。

向量化的处理方式是现代计算机的一个特点，无论硬件还是软件上，都提供了支持。

硬件上的支持：Intel’s MMX, SSE, and ARM’s NEON instruction sets support such vectorized loops.（摘自wikipedia：http://en.wikipedia.org/wiki/Vectorization_(computer_science)）

软件上的支持：著名的线性代数运算包blas对矩阵运算的自动并行化。例如，在R里执行如下的简单命令，在不需要作任何额外工作的情况下，就可以自动地实现并行化（前提是你的机器安装了blas）。

1
2

M = matrix(5000*5000, 5000, 5000)  ## 生成一个5000*5000的矩阵
S = M %*% M    ## 作矩阵乘法，在多核并安装了blas的机器里，会自动并行

应用向量化的思维方式去解决问题：

当我们习惯了用C语言的单步循环的思想来思考问题的时候，要把思维切换到向量化的方式是件比较困难的事件，以下显示一个例子，看用向量化的思想是怎么去解决问题，这样的描述又是多么美。

任务：对于稀疏矩阵M，让其第i(i＝1…m)行的各非零元素减去某个值w[i]。

正常的循环式思维的解决方案比较直观，作两层循环，让第i行非零元素减去w[i]即可，但这样的操作在脚本语言特别是像R这样的科学计算语言里执行效率不高，而当你面对的是一个稀疏矩阵时，问题又会变得复杂。

向量化的解决方案如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

## 利用向量化的计算范式设计算法
## 实例：稀疏矩阵不同行的各元素减去某个值
 
library(Matrix)
generate <- function(nr, nc, sr){
  L = nr*nc
  M = Matrix(data=0, nrow=nr, ncol=nc)
  idx = sample(1:L, as.integer(sr*L))
  val = rnorm(sr*L, mean=1)
  M[idx] = val
  return(M)
}
 
M = generate(10000, 5000, 0.1)  ## 生成稀疏矩阵
V = rnorm(10000)  ## 矩阵各行减去的值
N = M
N@x = 1  ## 复制并设置矩阵
A = Diagonal(nrow(N), V) %*% N  ## 把V的值放到新矩阵的各行非零单元中
M@x = M@x - A@x  # M = M-A

除去数据准备阶段，真正用于解决问题的只有第16到19行的4行命令，十分短小精悍，但需要一定的线性代数知识去理解这个过程。

从向量化到并行计算：

除了计算机硬件与科学计算包的支持外，向量化计算还是并行计算的天然先驱，如果你向量化后的算法效率仍然不佳，可以考虑把它并行化。

R那浩如烟海的第三方包里提供了大量支持并行计算的包，这里选择的是snowfall这个包，它可以简单地构建一个计算集群，把计算并行分布到集群的不同节点进行计算。以下通过一个例子比较循环，向量化以及并行三种方式在效率上的差异。

实战案例：for循环，apply，snowfall（并行）的比较。

任务：对一个矩阵每列排序并取回前10个。

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

library(snowfall)
 
mysort <- function(x){
#    replicate(2, sort(x))
    return(sort(x)[1:10])
}
 
do_for <- function(M){
    ans = matrix(0, 10, ncol(M))
    for(i in 1:ncol(M)){
        ans[,i] = sort(M)[1:10]
    }
}
 
do_apply <- function(M){
    return(apply(M, 2, mysort))
}
 
do_snowfall <- function(M, ncl){
    sfInit(parallel=TRUE, cpus=ncl)  ##初始化集群环境
    ans <- sfApply(M, 2, mysort)  ##把任务分发到各个slave
    sfStop()  ##关闭集群
    return(ans)
}
 
#M = matrix(rnorm(10000000), 100, 100000)
#system.time(ans <- do_for(M))
#system.time(ans <- do_apply(M))
#system.time(ans <- do_snowfall(M, 4))  ##用4个核并行

在我的服务器里运行的结果是这样的，do_for循环基本不可用，do_apply可以在一分钟内运算完毕，而do_snowfall的时间为do_apply的一半。当取消第4行的注释，即增加各个子任务的计算负荷后，do_apply与do_snowfall的计算时间都增加，而do_snowfall的时间变为do_apply的三分之一。

结论1：向量化的计算方式比起传统的循环计算有极大的性能优势。

结论2：由于并行的过程为：master把任务分解，分发到多个slave进行运算，slave返回结果到master。所以多核计算并不一定比最优的单核计算快，要看性能的瓶颈在slave还是在master上。

结论3：适合并行/分布计算的情景主要有两种：一是各slave的计算耗费为瓶颈，并行到多核能减少计算时间，越是slave耗时型的计算并行收益越大；二是一台机器的内存不足以进行整体的计算，分布到多台机器计算能把内存占用分开，这种情况下即使分布计算比单机慢也是可以接受的。

与map-reduce扯扯关系：

Map-reduce是google三驾马车之一，试图把所有计算都转换为map和reduce两种基本的运算方式，而达到良好的可并行性，从而实现计算上的可扩展性。虽然并不是每个公司都能实现像map-reduce那样的巨型计算框架，但是如果你能熟练地使用向量化的思想设计你的算法，那其实就是一个map-reduce的超集。向量化计算这个编程范式比map-reduce要复杂，但本质上都是使用了两种基本的运算，如果你把apply想像为map，把矩阵或矢量的乘运算想像为reduce，它们就是一体的。你总可以把向量化的计算通过apply和矩阵运算来实现，如果有一天你拥有了一个类map-reduce的计算框架，根据它们的对应关系，你的算法迁移就会非常地方便，甚至，你可以在你已经实现的向量化算法集当中抽取出一些共性来搭建自己的map-reduce系统。

别忘了一切优化的终极准则：过早的优化是魔鬼。

并行计算或分布式计算是为数据量与计算量日益膨胀的情况下所必须考虑的，但它的逻辑结构与维护成本也要高得多，如果你的应用场景还没达到需要进行并行或者分布计算的程度，没必要淌这趟混水。但是，向量化计算的编程范式却是很重要的，一个用向量化计算方式实现的算法，不但在当时获得了效益，而且其可扩展性也必定是强悍的，因为正如前文所述，向量化计算就是并行计算的天然先驱。除此之外，经常用向量化的方式来思考你的问题，你一定能感受到一种整体之美。

后记：

以上内容是从我参加第三届中国R语言会议的讲稿中整理出来的，以讲向量化的思想为主，R语言为辅，当时的PPT可以从这里下载（PPT上内容比较少，可以下载下来看备注部分）。我一直认为，通过这样的讲演，除了布道，更重要的是能够对自己一段时间以来的思考与工作做一个总结，其间无论是从自身整理还是从他人的提问，自己从中都能有所收获。根据当场及事后一些听众的疑问，我后续整理了下面一些内容。

如果问题很难向量化怎么办：并不是说R中不能使用for式的循环，事实上只有一层的for循环，并且对性能没有太大的影响的话，是不需要作过多的优化的。如果是两层for循环，可以考虑用apply去掉一层，如果有三层的循环，那很有可能就是你的算法有问题了。如果你的问题确实很难转换成向量化的形式，又或者说你已经懒于这样去做，那么不妨把最难于转换的一部分用C来实现，然后给R做一个接口，剩下的再考虑向量化会简单很多。我提倡的是结合R自有函数与向量化思想进行编程，因为R的自有函数大多都是用底层语言实现好的，效率有保证。

有人疑问是不是所有算法都可以并行：非基本的算法或多或少都是可以并行的，像kmeans，虽然从整体来说是一个迭代依赖的非并行式算法，但在每一步的迭代中却是可以实施并行的。应用map-reduce框架基于这么一个假设，你的问题，总可以找到一种算法来实现，这种算法可以或大部分可以用map加reduce基本操作来实现。如果向量化能跟map-reduce对应起来，那么也可以认为，你的问题，总能找到一种可以用向量化思想表示的算法来完全或部分解决，如果有些部分很难向量化，参考上一个问题。

我选择的算法本身就很复杂：有这么一个说法，如果你的数据不足够多，信息不够完整，你会想出很多奇思怪招来从中获得结论，实际上很多科学上精巧而复杂的算法是基于这种情况而设计出来的。而在实际中如果你的数据急剧膨胀，信息与噪声都足够地多，你的问题焦点就变为如何用一个或一些简单有效算法或算法的组合来提取有价值的信息。如果你在小数据集的时候习惯了采用复杂的算法，当你的应用场景转变了，数据规模变了，却仍旧沿用原有的思维方式，是注定要尴尬的（仅针对工程应用环境，而非科学应用环境）。

本文链接：http://www.wentrue.net/blog/?p=945

最初接触到足球时，94年世界杯已经远去，所以我只能一次次地从队友的口中重温巴乔那落寞的背影。那时，即使失败那也是一种悲情的浪漫。

98法兰西之夏，战火初燃，正好赶上我在积极地备战中考，绿茵场上的球星们在努力地追逐着大力神杯，而我也正希冀着找回自己三年前失落的那所重点中学。天气很炎热，每天晚自习后，我都会走到学校门口的小店里，要一碗豆腐脑，顺便瞟几眼屏幕里已经稍有认识的球员，然后匆匆赶回宿舍休息或是夜战。那时，我的眼中其实只有中国国家队，97年的金州兵败，是我第一次体验了作为中国球迷的苦涩。记得那时我完整地看完的一场小组赛是阿根廷对日本，因为那时我认识巴蒂、奥特加，还有传球失误导致失利的日本天才中场中田英寿。当然，那时我最喜欢的是无所不能的外星人罗纳尔多（那时我们叫他朗拿度）。及至中考结束，在家里，我第一次熬夜看球，看着巴西杀入决赛，然后金杯旁落，我第一次体会到在足球场上，球星不并等同于胜利。我也开始明白，作为一个球迷或是一名球员，失意是一件随时都要准备面对的事情。那时，我正逐步走出自己的童话王国。

02日韩世界杯时，我已经在读大二。由于宿舍没有电视，当时网络看直播也不成熟，我和同学相约跑到隔壁邮科院的电影院去看球，因为电影院不放电影了，专门摆了一个电子屏幕放球赛。于是我第一次有了全场人一起欢呼一起叹息的团体看球体验。其实只去过一次，因为第一次聚众大屏幕看球的体验就是眼睁睁地看着阿根廷负于英格兰，后来小组出局。随后的我，就像生活中大多数时候的状态一样，再也找不到做一件事情的理由与目的。所以，那届世界杯，除了中国队的比赛，我看过的赛事只有寥寥几场，也许是看似空闲的大学光阴，被我挥霍到网络与游戏世界中去了吧，抛却了曾经的钟爱与梦想，去追逐一些虚拟的快乐。那段时间，是我人生的第一段黑暗时期。但是，德国人的顽强，让我看到了一丝力量。这个衰落中的帝国，这个近年来已经越来越难找到闪耀巨星的日耳曼战车，每一次的关键赛事，却总能勇往直前。

06年的德国世界杯是迄今为止我看过最完整的一届世界杯，基本上绝大部分的比赛我都通过直播或录像看过。当时正是研一的尾声，那通常也是一个中国学生在学校系统学习的最后一段时间。基本上，我可以静下心来，做一些自己想做的事了。那一次我才了解到，原来世界杯最精彩的不是决赛，而是小组赛。当全世界球迷的目光都聚焦在那个沸腾的国家，当所有人都还有希望还有期待的时候，才是最美好的时刻。当几家欢笑几家愁一幕幕上映，当一个个希望化为一张张回国的机票，这个全球球迷的嘉年华就只退化为那座金杯。虽然最后捧杯的是铜墙铁壁的意大利，但我同样记住了华丽的阿根廷，激情四射的克林斯曼，逆境还魂的法国与悲情谢幕的齐祖，甚至灵魂附体的黄健翔，还有很多很多，有了这些，对于我们大部分人来说，远比那座金杯更来得珍贵。那时，我开始明白世界杯的真正意义。也是那年，虽然仍旧磕磕碰碰，我却意气风发。一个事物之所以精彩，是因为它能给所有人以希望与快乐，而非它的结局。

2010，又会上演什么样的故事？

———————————————————————————————————-

一直觉得，踢球的孩子更容易长大，因为他们很快地就能找到自己的队友，并且一开始就得在场上学习寻找自己合适的位置，继而扮演好自己的角色。无论成功或者失败的体验，他们都能首先在球场上获得。

球场是人生舞台的一个剪影，曾经年轻的我们都幻想过有朝一日自己能成为顶级俱乐部的球星，国家队的栋梁之才，尽情的在广阔的绿茵场上奔跑，胜利了撒足狂奔，失败了悲情落泪。到今天，我们默默地奔跑在离自己距离最近的一支小球队中，追逐着最简单的快乐，即使有时泛起那么点功利心，跟以前相比也渺小得那么地不值一提。

———————————————————————————————————-

到目前为止的这个夏天已经带我们太多的悲喜剧，赛季初挥金如土的银河战舰在四大皆空的同时，凭借着自身输出的斯内德和罗本，把国米和拜仁送进了欧冠决赛，行为艺术家伊布拉同学，赛季前为了欧冠离开国米，然后国米收获了40多年来从未染指的欧冠奖杯。我们没有理由不期待接下来更多的精彩。

最后，不能免俗，谈谈我看好的几支球队吧。

巴西：一旦巴西务实起来，你没有理由不看好他。

西班牙：西班牙依然看上去很美，如果他夺冠，估计是大多数“真球迷”所乐意看到的。

英格兰：英格兰的阵容其实很大牌的，还有一个那么大牌的教练。足球鼻祖，那么多年了，再拿一回大家估计也可以接受。

德国：哎，悲情的巴拉克，强悍的德国战车。

阿根廷：支持阿根廷，马拉多纳你滚蛋吧！

很累，他紧紧地依着自己的登山索，他知道一刻都不能放手，但真的很想睡上一觉。

嗡… … 耳边一直传来那种不知道来自哪个世界的杂音。

还好，他还有力气，接着走这段必然会有尽头的路。

轮回，是历史的一部分，也是生命的基本。你也许早已熟知，但它仍然会很残酷。