关联规则：R 与 SAS 的比较

啤酒和尿布的故事是关联分析方法最经典的案例，而用于关联分析的 Apriori 算法更是十大数据挖掘算法之一（http://www.cs.uvm.edu/~icdm/algorithms/index.shtml，这个排名虽然是几年前的调查结果，但是其重要性仍可见一斑）。本文以《R and Data Mining》书中使用的泰坦尼克号人员的生存数据为例，介绍如何使用 R 和 SAS 的 Apriori 算法进行关联分析，比较两者的建模结果并对结果中存在的差异进行解释分析。

一、关联分析

网上有很多资料介绍关联分析算法，本文就不再赘述。我自己看的是《Introduction to Data Mining》(有对应的中文版，人民邮电出版社的《数据挖掘导论》)，愿意看英文的同学可以访问：http://www-users.cs.umn.edu/~kumar/dmbook/ch6.pdf。网上其他的资料我也大致翻过，对比之后感觉这本书是一本相当不错的教材，算法方面介绍地比较全面且有一定深度。我本人不建议大家去看那些非专业人士总结的关联分析算法介绍，虽然浅显易懂，但是内容片面，容易误导初学者，错把树木当成了森林。

对于关联分析在行业应用中的经验分享、初学者的误区和最佳实践方面的资料很少，唯一能找到的一本好书是清华大学出版社的《啤酒与尿布》，主要介绍购物篮分析在零售行业的应用。我始终认为分析师除了算法和软件，还需要了解行业背景，不然挖出的只是模式，而不是切实可行并且能带来商业价值的模式，甚至还有可能是错误的模式。 

二、软件

我只用过 R 和 SAS，其他的软件没碰过，所以只能对这两个软件进行比较。

<td>
  <strong>R/ARULES</strong>
</td>

<td>
  <strong>SAS/EM</strong>
</td>

<td>
  Yes
</td>

<td>
  Yes
</td>

<td>
  Yes
</td>

<td>
  No
</td>

<td>
  No
</td>

<td>
  No
</td>

据网友说 Excel 也能做关联分析，但是因为其对数据进行了抽样，所以每次运行的结果都不一样。SPSS 的 Modeler 不知道怎么样，有用过的同学请分享一下经验，最好使用泰坦尼克号的数据进行分析，这样可以比较一下各软件的结果是否相同。

三、R 的代码和结果

R 的代码主要来自《R and Data Mining》，我只加了下载数据的代码和对代码的中文说明。

1）下载泰坦尼克数据

setInternet2(TRUE)

con <- url(“http://www.rdatamining.com/data/titanic.raw.rdata“)

load(con)

close(con) # url() always opens the connection

str(titanic.raw)

2）关联分析

library(arules)

# find association rules with default settings

rules <- apriori(titanic.raw)

inspect(rules)

3）只保留结果中包含生存变量的关联规则

# rules with rhs containing “Survived” only

rules <- apriori(titanic.raw, parameter = list(minlen=2, supp=0.005, conf=0.8), appearance = list(rhs=c(“Survived=No”, “Survived=Yes”), default=“lhs”),control = list(verbose=F))

rules.sorted <- sort(rules, by=“lift”)

inspect(rules.sorted)

R 总共生成了 12 条跟人员生存相关的规则：

lhs       rhs      support      confidence      lift

1 {Class=2nd, Age=Child}                         => {Survived=Yes}

0.010904134 1.0000000 3.095640

2 {Class=2nd, Sex=Female, Age=Child}  => {Survived=Yes}

0.005906406 1.0000000 3.095640

3 {Class=1st, Sex=Female}                      => {Survived=Yes}

0.064061790 0.9724138 3.010243

4 {Class=1st, Sex=Female, Age=Adult}    => {Survived=Yes}

0.063607451 0.9722222 3.009650

5 {Class=2nd, Sex=Male, Age=Adult}        => {Survived=No}

0.069968196 0.9166667 1.354083

6 {Class=2nd, Sex=Female}                      => {Survived=Yes}

0.042253521 0.8773585 2.715986

7 {Class=Crew, Sex=Female}                   => {Survived=Yes}

0.009086779 0.8695652 2.691861

8 {Class=Crew, Sex=Female, Age=Adult} => {Survived=Yes}

0.009086779 0.8695652 2.691861

9 {Class=2nd, Sex=Male}                           => {Survived=No}

0.069968196 0.8603352 1.270871

10 {Class=2nd, Sex=Female, Age=Adult}  => {Survived=Yes}

0.036347115 0.8602151 2.662916

11 {Class=3rd, Sex=Male, Age=Adult}       => {Survived=No}

0.175829169 0.8376623 1.237379

12 {Class=3rd, Sex=Male}                          => {Survived=No}

0.191731031 0.8274510 1.222295

4）去除冗余的规则

# find redundant rules

subset.matrix <- is.subset(rules.sorted, rules.sorted)

subset.matrix[lower.tri(subset.matrix, diag=T)] <- NA

redundant <- colSums(subset.matrix, na.rm=T) >= 1

which(redundant)

# remove redundant rules

rules.pruned <- rules.sorted[!redundant]

inspect(rules.pruned)

去除冗余的规则后剩下 8 条规则：

lhs       rhs      support      confidence      lift

1 {Class=2nd, Age=Child}                   => {Survived=Yes}

0.010904134  1.0000000 3.095640

2 {Class=1st, Sex=Female}                => {Survived=Yes}
0.064061790  0.9724138 3.010243
3 {Class=2nd, Sex=Female}               => {Survived=Yes}
0.042253521  0.8773585 2.715986
4 {Class=Crew, Sex=Female}            => {Survived=Yes}

0.009086779  0.8695652 2.691861

5 {Class=2nd, Sex=Male, Age=Adult} => {Survived=No}

0.069968196  0.9166667 1.354083

6 {Class=2nd, Sex=Male}                   => {Survived=No}

0.069968196  0.8603352 1.270871

7 {Class=3rd, Sex=Male, Age=Adult}  => {Survived=No}

0.175829169  0.8376623 1.237379

8 {Class=3rd, Sex=Male}                    => {Survived=No}

0.191731031  0.8274510 1.222295

5）结果的解释

对于结果的解释，一定要慎重，千万不要盲目下结论。从下面的四条规则看，好像确实像电影中描述的那样：妇女和儿童优先。

1 {Class=2nd, Age=Child}              => {Survived=Yes} 0.010904134  1.0000000 3.095640

2 {Class=1st, Sex=Female}           => {Survived=Yes} 0.064061790  0.9724138 3.010243

3 {Class=2nd, Sex=Female}          => {Survived=Yes} 0.042253521  0.8773585 2.715986

4 {Class=Crew, Sex=Female}       => {Survived=Yes} 0.009086779  0.8695652 2.691861

如果我们减小最小支持率和置信度的阈值，则能看到更多的真相。

rules <- apriori(titanic.raw, parameter = list(minlen=3, supp=0.002, conf=0.2), appearance = list(rhs=c(“Survived=Yes”), lhs=c(“Class=1st”, “Class=2nd”, “Class=3rd”, “Age=Child”, “Age=Adult”), default=“none”), control = list(verbose=F))

rules.sorted <- sort(rules, by=“confidence”)

inspect(rules.sorted)

lhs                        rhs           support     confidence lift

1 {Class=2nd, Age=Child} => {Survived=Yes} 0.010904134 1.0000000 3.0956399

2 {Class=1st, Age=Child} => {Survived=Yes} 0.002726034 1.0000000 3.0956399

3 {Class=1st, Age=Adult} => {Survived=Yes} 0.089504771 0.6175549 1.9117275
4 {Class=2nd, Age=Adult} => {Survived=Yes} 0.042707860 0.3601533 1.1149048

5 {Class=3rd, Age=Child} => {Survived=Yes} 0.012267151 0.3417722 1.0580035
6 {Class=3rd, Age=Adult} => {Survived=Yes} 0.068605179 0.2408293 0.7455209

从规则 3 和规则 5 以及之前的规则 2 和 3 可以看出泰坦尼克号获得优先权的主要是头等舱、二等舱的妇孺。

据统计，头等舱男乘客的生还率比三等舱中儿童的生还率还稍高一点。美国新泽西州州立大学教授、著名社会学家戴维 · 波普诺研究后毫不客气地修改了曾使英国人颇感 “安慰” 的“社会规范”(妇女和儿童优先)：“在泰坦尼克号上实践的社会规范这样表述可能更准确一些：‘头等舱和二等舱的妇女和儿童优先’。”

这些是关于泰坦尼克号生存数据分析的资料：

泰坦尼克号逃生真相：“妇女儿童优先” 只是个传说

历史没有那么温暖

6）可视化

# visualize rules

library(arulesViz)

plot(rules)

plot(rules, method=“graph”, control=list(type=“items”))

plot(rules, method=“paracoord”, control=list(reorder=TRUE))

对于不熟悉 R 的 SAS 用户，可以阅读以下资料学习 R 以及 ARULES 包：

http://cran.r-project.org/web/packages/arules/vignettes/arules.pdf

https://science.nature.nps.gov/im/datamgmt/statistics/R/documents/R_for_SAS_SPSS_users.pdf

四、SAS 代码和结果

1）下载泰坦尼克数据

proc iml;
submit /R;
setInternet2(TRUE)

con <- url(http://www.rdatamining.com/data/titanic.raw.rdata)
load(con)
close(con) # url() always opens the connection

endsubmit;

call ImportDataSetFromR(“Work.titanic”, “titanic.raw”);
run;quit;

2）将数据转换成 SAS/EM 要求的格式

data items2;
set titanic;
length tid 8;
length item $8;
tid = _n_;

item = class;

output;
item = sex;

output;
item = age;

output;
item = survived;

output;
keep tid item;
run;

3）关联分析

proc dmdb data=items2 dmdbcat=dbcat;
class tid item;
run; quit;

proc assoc data=items2 dmdbcat=dbcat pctsup=0.5 out=frequentItems;
id tid;
target item;

run;

proc rulegen in=frequentItems dmdbcat=dbcat out=rules minconf=80;
run ;

proc sort data=rules;
by descending conf;
run ;

4） 只保留结果中包含生存变量的关联规则

data surviverules;
set rules(where=(set_size>1 and (_rhand=‘Yes’ or _rhand=‘No’)));
run;

proc print data=surviverules;
var conf support lift rule ;
run ;

SAS 结果:

<th scope="col">
  CONF
</th>

<th scope="col">
  SUPPORT
</th>

<th scope="col">
  LIFT
</th>

<th scope="col">
  RULE
</th>

<td>
  100.00
</td>

<td>
  1.09
</td>

<td>
  3.10
</td>

<td>
  2nd & Child ==> Yes
</td>

<td>
  100.00
</td>

<td>
  0.59
</td>

<td>
  3.10
</td>

<td>
  2nd & Child & Female ==> Yes
</td>

<td>
  100.00
</td>

<td>
  0.50
</td>

<td>
  3.10
</td>

<td>
  2nd & Child & Male ==> Yes
</td>

<td>
  97.24
</td>

<td>
  6.41
</td>

<td>
  3.01
</td>

<td>
  1st & Female ==> Yes
</td>

<td>
  97.22
</td>

<td>
  6.36
</td>

<td>
  3.01
</td>

<td>
  1st & Adult & Female ==> Yes
</td>

<td>
  91.67
</td>

<td>
  7.00
</td>

<td>
  1.35
</td>

<td>
  2nd & Adult & Male ==> No
</td>

<td>
  87.74
</td>

<td>
  4.23
</td>

<td>
  2.72
</td>

<td>
  2nd & Female ==> Yes
</td>

<td>
  86.96
</td>

<td>
  0.91
</td>

<td>
  2.69
</td>

<td>
  Crew & Female ==> Yes
</td>

<td>
  86.96
</td>

<td>
  0.91
</td>

<td>
  2.69
</td>

<td>
  Adult & Crew & Female ==> Yes
</td>

<td>
  86.03
</td>

<td>
  7.00
</td>

<td>
  1.27
</td>

<td>
  2nd & Male ==> No
</td>

<td>
  86.02
</td>

<td>
  3.63
</td>

<td>
  2.66
</td>

<td>
  2nd & Adult & Female ==> Yes
</td>

<td>
  83.77
</td>

<td>
  17.58
</td>

<td>
  1.24
</td>

<td>
  3rd & Adult & Male ==> No
</td>

<td>
  82.75
</td>

<td>
  19.17
</td>

<td>
  1.22
</td>

<td>
  3rd & Male ==> No
</td>

有关 SAS/EM 关联分析的公开资料很少，产品的在线帮助文档大概从 4.3 以后的版本就设置了访问权限，只有 SAS/EM 的用户才能阅读，新版本的功能和界面跟 4.3 版本有很大差别。这里只能给大家一些 4.3 的帮助文档，主要是上面代码中用到的几个过程步：http://support.sas.com/documentation/onlinedoc/miner/em43/dmdb.pdf

http://support.sas.com/documentation/onlinedoc/miner/em43/assoc.pdf

http://support.sas.com/documentation/onlinedoc/miner/em43/sequence.pdf

http://support.sas.com/documentation/onlinedoc/miner/em43/rulegen.pdf

mbscore(购物篮数据的预测，是 EM 6.1/SAS 9.2 时新引入的过程步，支持层次关联)

五、结果比较

从上面的结果看，R 生成了 12 条规则，而 SAS 生成了 13 条规则，对比每条规则后，发现 SAS 的第 3 条规则在 R 中没有。

<td>
  100.00
</td>

<td>
  0.50
</td>

<td>
  3.10
</td>

<td>
  2nd & Child & Male ==> Yes
</td>

我猜测原因是两个软件对最小支持度的处理不太一样，SAS 可能是对最小支持度百分比乘以总记录条数后取整了。此处，泰坦尼克数据总共有 2201 条记录，最小支持度百分比为 0.5%，两者相乘积为 11.005，而 2nd & Child & Male ==> Yes 这条规则总共出现过 11 次，如果严格按照实数大小比较，不应该出现在最后的结果中，但是如果按照整数部分比较，则结果正确。打算将 SAS 模型切换到 R 或者将 R 模型切换到 SAS 的同学要注意这个差异，结果有时不完全一样！

data min_support;
set frequentItems;
if count=int(2201*0.005);
run;

proc print data=min_support;
run;quit; 

<th scope="col">
  SET_SIZE
</th>

<th scope="col">
  COUNT
</th>

<th scope="col">
  ITEM1
</th>

<th scope="col">
  ITEM2
</th>

<th scope="col">
  ITEM3
</th>

<th scope="col">
  ITEM4
</th>

<th scope="col">
  ITEM5
</th>

<th scope="col">
  ITEM6
</th>

<td>
  3
</td>

<td>
  11
</td>

<td>
  2nd
</td>

<td>
  Child
</td>

<td>
  Male
</td>

<td>
   
</td>

<td>
   
</td>

<td>
   
</td>

<td>
  4
</td>

<td>
  11
</td>

<td>
  2nd
</td>

<td>
  Child
</td>

<td>
  Male
</td>

<td>
  Yes
</td>

<td>
   
</td>

<td>
   
</td>

相比 SAS，R 关联分析中比较吸引人的功能就是从规则集中去除冗余的规则，这一功能 SAS 里面好像没有（我没找到）。SAS 用户如果想要使用 R 的这个功能，我找到的唯一办法就是将 SAS 的关联规则导出成 PMML 文件，然后再将 PMML 文件导入 R 生成对应的 Rule 对象，但是这个方法因为我的环境有点问题，所以我自己没试。

有兴趣的同学，可以看看下面的资料：

1）如何将 PMML 文件导入 R 生成 Rule 对象？

2）如何在 SAS EMM 中使用 PMML？

附：PMML 技术的未来

对于模型的部署和使用，尤其是跨软件、平台的使用场景下或者对于大数据的分析，PMML 是一个可行的解决方案，有一些厂商已经在自己的产品中通过 PMML 这种方式来实现对大数据的分析预测。

Zementis：

Deploying Predictive Analytics with PMML, R evolution R, and ADAPA

PMML: Accelerating the Time to Value for Predictive Analytics in the Big Data Era

IBM：Database Mining Guide