用R语言进行关联分析

关联是两个或多个变量取值之间存在的一类重要的可被发现的某种规律性。关联分析目的是寻找给定数据记录集中数据项之间隐藏的关联关系，描述数据之间的密切度。

几个基本概念

1. 项集

这是一个集合的概念，在一篮子商品中的一件消费品即为一项（Item），则若干项的集合为项集，如{啤酒，尿布}构成一个二元项集。

2. 关联规则

一般记为的形式，X为先决条件，Y为相应的关联结果，用于表示数据内隐含的关联性。如：，表示购买了尿布的消费者往往也会购买啤酒。

关联性强度如何，由三个概念——支持度、置信度、提升度来控制和评价。

例：有10000个消费者购买了商品，其中购买尿布1000个，购买啤酒2000个，购买面包500个，同时购买尿布和面包800个，同时购买尿布和面包100个。

3. 支持度（Support）

支持度是指在所有项集中{X, Y}出现的可能性，即项集中同时含有X和Y的概率：

该指标作为建立强关联规则的第一个门槛，衡量了所考察关联规则在“量”上的多少。通过设定最小阈值（minsup），剔除“出镜率”较低的无意义规则，保留出现较为频繁的项集所隐含的规则。

设定最小阈值为5%，由于{尿布，啤酒}的支持度为800/10000=8%，满足基本输了要求，成为频繁项集，保留规则；而{尿布，面包}的支持度为100/10000=1%，被剔除。

4. 置信度（Confidence）

置信度表示在先决条件X发生的条件下，关联结果Y发生的概率：

这是生成强关联规则的第二个门槛，衡量了所考察的关联规则在“质”上的可靠性。相似的，我们需要对置信度设定最小阈值（mincon）来实现进一步筛选。

具体的，当设定置信度的最小阈值为70%时，置信度为800/1000=80%，而的置信度为800/2000=40%，被剔除。

5. 提升度（lift）

提升度表示在含有X的条件下同时含有Y的可能性与没有X这个条件下项集中含有Y的可能性之比：

该指标与置信度同样衡量规则的可靠性，可以看作是置信度的一种互补指标。

R中Apriori算法

算法步骤：

1. 选出满足支持度最小阈值的所有项集，即频繁项集；

2. 从频繁项集中找出满足最小置信度的所有规则。

>library(arules) #加载arules包

>click_detail =read.transactions("click_detail.txt",format="basket",sep=",",cols=c(1)) #读取txt文档（文档编码为ANSI）

>rules <- apriori(click_detail, parameter =list(supp=0.01,conf=0.5,target="rules")) #调用apriori算法

>rules

set of419 rules

>inspect(rules[1:10]) #查看前十条规则

解释

1)library(arules)：加载程序包arules，当然如果你前面没有下载过这个包，就要先install.packages(arules)

2)click_detail =read.transactions("click_detail.txt",format="basket",sep=",",cols=c(1))：读入数据

read.transactions(file, format =c("basket", "single"), sep = NULL,

cols = NULL, rm.duplicates =FALSE, encoding = "unknown")

file:文件名，对应click_detail中的“click_detail.txt”

format:文件格式，可以有两种，分别为“basket”,“single”，click_detail.txt中用的是basket。

basket: basket就是篮子，一个顾客买的东西都放到同一个篮子，所有顾客的transactions就是一个个篮子的组合结果。如下形式，每条交易都是独立的。

文件形式：

item1,item2

item1

item2,item3

读入后：

items

1 {item1,

item2}

2 {item1}

3 {item2,

item3}

single: single的意思，顾名思义，就是单独的交易，简单说，交易记录为：顾客1买了产品1， 顾客1买了产品2，顾客2买了产品3……(产品1，产品2，产品3中可以是单个产品，也可以是多个产品)，如下形式：

trans1 item1

trans2 item1

trans2 item2

读入后：

items transactionID

1 {item1}trans1

2 {item1,

item2}trans2

sep:文件中数据是怎么被分隔的，默认为空格，click_detail里面用逗号分隔

cols:对basket, col=1,表示第一列是数据的transaction ids(交易号)，如果col=NULL，则表示数据里面没有交易号这一列；对single，col=c(1,2)表示第一列是transaction ids，第二列是item ids

rm.duplicates:是否移除重复项，默认为FALSE

encoding:写到这里研究了encoding是什么意思，发现前面txt可以不是”ANSI”类型，如果TXT是“UTF-8”，写encoding=”UTF-8”，就OK了.

3)rules <- apriori(click_detail,parameter = list(supp=0.01,conf=0.5,target="rules"))：apriori函数

apriori(data, parameter = NULL, appearance = NULL, control = NULL)

data:数据

parameter：设置参数，默认情况下parameter=list(supp=0.1,conf=0.8,maxlen=10,minlen=1,target=”rules”)

supp:支持度（support）

conf:置信度（confidence）

maxlen,minlen:每个项集所含项数的最大最小值

target:“rules”或“frequent itemsets”（输出关联规则/频繁项集）

apperence:对先决条件X（lhs），关联结果Y（rhs）中具体包含哪些项进行限制，如：设置lhs=beer，将仅输出lhs含有beer这一项的关联规则。默认情况下，所有项都将无限制出现。

control：控制函数性能，如可以设定对项集进行升序sort=1或降序sort=-1排序，是否向使用者报告进程（verbose=F/T）

补充

通过支持度控制：rules.sorted_sup = sort(rules, by=”support”)

通过置信度控制：rules.sorted_con = sort(rules, by=”confidence”)

通过提升度控制：rules.sorted_lift = sort(rules, by=”lift”)

Apriori算法

两步法：

1. 频繁项集的产生：找出所有满足最小支持度阈值的项集，称为频繁项集；

2. 规则的产生：对于每一个频繁项集l，找出其中所有的非空子集；然后，对于每一个这样的子集a，如果support(l)与support(a)的比值大于最小可信度，则存在规则a==>(l-a)。

频繁项集产生所需要的计算开销远大于规则产生所需的计算开销

频繁项集的产生

几个概念：

1， 一个包含K个项的数据集，可能产生2^k个候选集

2，先验原理：如果一个项集是频繁的，则它的所有子集也是频繁的（理解了频繁项集的意义，这句话很容易理解的）；相反，如果一个项集是非频繁的，则它所有子集也一定是非频繁的。

3基于支持度（SUPPORT）度量的一个关键性质：一个项集的支持度不会超过它的子集的支持度（很好理解，支持度是共同发生的概率，假设项集{A,B,C}，{A,B}是它的一个自己，A,B,C同时发生的概率肯定不会超过A,B同时发生的概率）。

上面这条规则就是Apriori中使用到的，如下图，当寻找频繁项集时，从上往下扫描，当遇到一个项集是非频繁项集（该项集支持度小于Minsup），那么它下面的项集肯定就是非频繁项集，这一部分就剪枝掉了。

一个例子（百度到的一个PPT上的）：

当我在理解频繁项集的意义时，在R上简单的复现了这个例子，这里采用了eclat算法，跟apriori应该差不多：

代码：

item <- list(

c("bread","milk"),

c("bread","diaper","beer","eggs"),

c("milk","diaper","beer","coke"),

c("bread","milk","diaper","beer"),

c("bread","milk","diaper","coke")

)

names(item) <- paste("tr",c(1:5),sep = "")

item

trans <- as(item,"transactions") #将List转为transactions型

rules = eclat(trans,parameter = list(supp = 0.6,

target ="frequent itemsets"),control = list(sort=1))

inspect(rules) #查看频繁项集

运行后结果：

>inspect(rules)

items support

1{beer,

diaper}0.6

2{diaper,

milk} 0.6

3{bread,

diaper}0.6

4{bread,

milk} 0.6

5{beer} 0.6

6{milk} 0.8

7{bread} 0.8

8{diaper} 0.8

以上就是该例子的所有频繁项集，然后我发现少了{bread,milk,diaper}这个项集，回到例子一看，这个项集实际上只出现了两次，所以是没有这个项集的。

规则的产生

每个频繁k项集能产生最多2k-2个关联规则

将项集Y划分成两个非空的子集X和Y-X，使得X ->Y-X满足置信度阈值

定理：如果规则X->Y-X不满足置信度阈值，则X’->Y-X’的规则一定也不满足置信度阈值，其中X’是X的子集

Apriori按下图进行逐层计算，当发现一个不满足置信度的项集后，该项集所有子集的规则都可以剪枝掉了。

本文引入了一个新的Stata命令xtpmg，用于估计包含大N大T的非平稳非平衡面板。基于非平稳面板文献的最新进展，xtpmg提供了三种替代估计量：1）传统的固定效应（FE）、2）Pesaran和Smith的平均组估计量MG（估计动态异质面板的长期关系）；3）Pesaran、Shin和Smith的混合平均组估计PMG（估计动态异质面板中的长期关系）。

摘要

一、引言

二、MG和PMG的估计原理

三、xtpmg命令在stata中的使用

    （一）基本句式

    （二）详细说明

四、实证操作

    （一）模型设定

    （二）PMG估计

    （三）MG估计

    （四）动态固定效应

五、结论

参考文献

近年来，与动态面板数据相关的文献开始关注截面观测数（N）和时间序列观测数（T）都较大的面板。海量数据的可获得性无疑是导致这种转变的关键因素。例如，一些跨国家或省份的数据集现在已经足够大，以至于可以分别估计每个国家（或省份）的参数。

大N大T动态面板的渐近性不同于传统的大N小T动态面板的渐近性。小T面板估计通常依赖于固定效应或随机效应估计，或固定效应估计和工具变量估计的组合，如Arellano和Bond（1991）的广义矩估计方法（GMM）。这些方法需要汇集（pool）各个组，并且只允许截距项在组间不同。然而，大N大T文献的一个中心发现是，斜率参数均匀性的假设通常是不合适的。Pesaran和Smith（1995年）、Im、Pesaran和Shin（2003年）、Pesaran、Shin和Smith（1997年、1999年）以及Phillips和Moon（2000年）提出了这一点。

随着大N大T动态面板固有的时间观测值的增加，非平稳性也受到关注。Pesaran，Shin和Smith（1997，1999）最近的论文提出了两种重要的新方法来估计非平稳动态面板：平均组（MG）和混合平均组（PMG）估计，其允许参数在不同组之间不同。MG估计值（参见Pesaran和Smith 1995）依赖于估计N个时间序列回归并平均系数，而PMG估计（参见Pesaran、Shin和Smith 1997、1999）依赖于系数的合并（pool）和平均。

在最近的实证研究中，MG和PMG估计量被应用于各种情形。例如，Freeman（2000）使用这些方法来评估1961-1995年间各州的酒精消费量。Martinez-Zarzoso和Bengochea-Morancho（2004）在1975-1998年间对22个经合组织成员国的环境库兹涅茨曲线进行了估算。Frank（2005）使用MG和PMG估计值来评估1945-2001年间美国各州收入不平等对经济增长的长期影响。

假设自回归分布滞后（ARDL）( )的动态面板模型具有如下形式：

其中横截面标记（组数）为i=1,2，…，N；时间序列标记（周期数）为t=1,2，…，T； 是解释变量，为k×1向量； 是k×1系数向量； 是标量； 是组特定效应。T必须足够大，以便模型可以分别估计每个组。模型也可以包括时间趋势和其他固定效应。

如果公式（1）中的变量是 I(1) 和协整的，那么误差项是所有 的 I(0) 过程。协整变量的一个主要特征是它们对长期均衡的任何偏离的响应。这一特征意味着在误差修正模型中，变量的短期动态变化受其偏离平衡的影响。因此，通常将公式（1）再次参数化为误差校正（Error Correction, EC）方程：

参数 是调整项的误差修正速度。如果 =0，则没有证据表明变量间存在长期关系。假设变量显示回归长期均衡的情况下，该参数预计将显著为负。特别重要的是向量 ，它包含变量之间的长期关系。

最近关于大N和大T的动态非平衡面板估计的文献提出了估计公式（2）的几种方法。在一种极端情况下，可以使用固定效应（FE）估计方法，其中每个组的时间序列数据被汇集（pool）在一起，并且只允许截距在组间不同。然而，如果斜率系数事实上不相同，那么FE方法会产生不一致且可能误导的结果。另一个极端情况下，Pesaran和Smith（1995）提出的MG估计可以分别为每组拟合，并可以计算出系数的简单算术平均值。使用这种估计方法，截距、斜率系数和误差方差都允许组间不同。

最近，Pesaran、Shin和Smith（1997、1999）提出了一个PMG估计法，它结合了合并和平均（pooling and averaging）。这种中间估计方法允许截距、短期系数和误差方差在组间不同（与MG估计法一样），但限制长期系数在组间相等（与FE估计法一样）。由于公式（2）的参数是非线性的，Pesaran、Shin和Smith（1999）提出了一种估计参数的极大似然法（Maximum likelihood method）。

将概率（likelihood）表示为每个横截面的概率（likelihood）的乘积并取对数，得到如下表达式：

xtpmg使用Stata强大的ML框架来实现PMG估计。具体来说，我们利用ml的hold选项，通过“反向替换”来最大化概率（likelihood）。从长期系数向量 的初始估计开始，可通过 对（ , ）的回归来估计短期系数和群组层面（group -specific）的调整项的速度。这些条件估计又被用来更新θ的估计。迭代过程一直持续直到收敛。

迭代条件似然最大化的参数估计与完全信息极大似然的参数估计是渐近一致的。但是，它们估计的协方差矩阵不同。然而，由于PMG参数的分布是已知的，我们可以获得所有估计参数的完全协方差矩阵。如Pesaran、Shin和Smith（1999）所示，协方差矩阵可以通过如下的逆矩阵求得：

MG参数只是单个系数的未加权平均值。例如，误差校正系数 的MG估计为

其方差为

其他短期系数的均值和方差也作了类似的估计。

lr（varlist） ——指定在计算长期协整向量时要包含的变量。

ec（string） ——用于指定新创建的误差更正的名称；默认值为_ _ec。

replace ——覆盖错误更正变量（如果存在）。

constraints( string ) ——指定要应用于模型的约束。此选项当前仅支持选项pmg一起使用。

noconstant ——不包含常数项。此选项不能与选项dfe一起使用。

cluster（ varname ） ——指定观察值在组间是独立的，但不一定在组内。varname指定每个观察所属的组，例如，对个体进行重复观察的数据中的cluster（personid）。cluster（）影响估计量的估计标准误差和方差-协方差矩阵（VCE），但不影响估计系数。

level(#) ——设置置信水平，默认是level(95)

technique(algorithm_spec) ——指定ml最大化方法。algorithm_spec是algorithm[#[algorithm[#]...]。algorithm可以是[nr|bfgs|dfp]。bhhh算法与xtpmg不兼容。technology（）只能与选项pmg一起使用。

difficult ——将使用不同的步进算法在非凹区域的概率（likelihood）。

full ——指定列出所有N个横截面回归结果。默认情况下，仅列出平均系数。

model ——是要拟合的估计方法，是以下类型之一：

    pmg 是默认值，并指定pmg估计。该模型限制长期系数向量在面板之间相等，同时允许组特定的短期和调整系数。

    mg 指定mg估计。该模型将参数拟合为N个个体组回归的平均值。

    dfe 指定了动态固定效应估计。

我们用24个OECD国家的年度总消费数据来说明xtpmg的使用。这些数据来自Pesaran、Shin和Smith（1997、1999），涵盖了1960-1993年。比利时1993年的年度观测数据不包括在估算样本中，比利时的估算期为1962-1992年，其他23个经合组织国家的估算期为1962-1993年。xtpmg要求在估计之前tsset设置数据。

假设长期消费函数为

其中，国家数i=1,2，…，N；周期数t=1,2，…，t；c是实际人均消费的对数；y是实际人均收入的对数；π是通货膨胀率。如果变量是I（1）和协整的，那么所有i的误差项都是I（0）过程。

公式（6）的ARDL（1,1,1）动态面板模型为

公式（7）的误差更正的再参数化公式是：

平差参数 和长期系数θ1i和θ2i的误差修正速度是最重要的。如果包含θ0i，则允许协整关系的非零均值。如果变量回归到长期均衡，人们会认为 为负。大多数总消费理论认为，长期收入弹性θ1i应等于1。通货膨胀效应θ2i通常被认为是负的。

首先，文章估计模型（8）的PMG估计量。在这种背景下，PMG模型考虑了异质的短期动态以及共同的长期收入和通货膨胀弹性。通常只有长期参数才有意义。pmg选项的默认结果包括长期参数估计和平均短期参数估计。

在结果中，估计的长期通胀弹性与预期一样显著为负。另外，估计的收入弹性也显著为正。理论上，收入弹性等于1。这个假设很容易验证：

相应的χ2值为121.2，可以拒绝单位收入弹性的原假设。

full选项估计并保存一个N+1多方程模型。第一个方程（按选项ec标注）表示标准化协整向量。剩下的N个方程列出了组特定的短期系数。

由于每组都有自己的估计方程，我们可以直观地预测变量。

类似地，交叉方程限制也很容易得到。

MG估计值是N个单独回归系数的未加权平均值。带有mg选项的xtpmg循环遍历样本中的所有面板，以估计（8）的参数

MG估计是一个两方程模型：标准化协整向量(EC)和短期动态系数（SR）。在比较PMG和MG估计量时，我们注意到，估计的长期收入和通货膨胀弹性在两个模型中都具有统计显著性和预期相同的符号。然而，PMG对通货膨胀弹性的估计在数量级上大于MG模型的估计（分别为-47和-35）。估计的长期收入弹性（分别为.90和.92）则正好相反。每个模型的调整估计速度意味着短期动态显著不同。

回想一下，PMG估计将所有面板的长期弹性限制为相等。当这些约束条件是真的时，这种跨国家的“加总”（pooling）产生了有效和一致的估计。然而，通常情况下，斜率相同的假设在经验上被否定。如果真实模型是异质的，则PMG估计是不一致的；MG估计在任何一种情况下都是一致的。这些模型的差异性检验是用常见的Hausman检验进行的。

计算的Hausman统计量为1.06，服从χ2（2）分布，应接受原假设。在这里我们得出结论，在原假设下（PMG和MG估计差异非系统性的），PMG估计是首选的有效估计。

动态FE估计与PMG估计一样，限制所有面板的协整向量系数相等。FE模型进一步限制了调整系数和短期系数的速度相等。带有dfe选项的xtpmg适合模型（8），同时允许面板层面的截距项。在计算标准误差时，允许使用cluster（）选项进行组内相关性计算。

动态FE模型的所有系数都产生了预期的符号，事实上，它们与PMG和MG估计值相似。正如Baltagi、Griffin和Xiong（2000）所讨论的，FE模型受到误差项和滞后因变量之间内生性的联立方程偏差的影响。豪斯曼检验（ Hausman test）可以很容易地用来衡量这种内生性的程度。

结果表明，在该样本数据下，联立方程偏差是最小的。在这个例子中，我们得出结论，FE模型优于MG模型。

本文介绍了Pesaran和Smith（1995）以及Pesaran，Shin和Smith（1997；1999）在估算具有大N大T的非平稳非平衡面板方面的最新进展。我们提供了一个新的Stata命令xtpmg，该命令估计了三种可供选择的模型：一个依赖于横截面合并(pooling)的传统动态FE估算，依赖于横截面平均值的MG估计和依赖于系数合并(pooling)和平均的PMG估计。

Blackburne III, E. F., &Frank, M. W. (2007). Estimation of nonstationary heterogeneous panels. The Stata Journal, 7(2), 197-208. 点击链接可在线阅读原文<https://maiimg.com/pdf/?e=agFqt4MbrUBcwm>

Arellano, M., and S. Bond. 1991. Some tests of specification for panel data: Monte Carlo evidence and an application to employment equations. Review of Economic Studies 58: 277–297.

Baltagi, B. H. 2001. Econometric Analysis of Panel Data. 2nd ed. New York: Wiley.

Baltagi, B. H., J. M. Griffin, and W. Xiong. 2000. To pool or not to pool: Homogeneous versus heterogeneous estimators applied to cigarette demand. Review of Economics and Statistics 82: 117–126.

Baum, C. F., M. E. Schaffer, and S. Stillman. 2003. Instrumental variables and GMM: Estimation and testing. Stata Journal 3: 1–31.

Frank, M. W. 2005. Income inequality and economic growth in the U.S.: A panel cointegration approach. Sam Houston State University Working Paper 05-03.

Freeman, D. G. 2000. Alternative panel estimates of alcohol demand, taxation, and the business cycle. Southern Economic Journal 67: 325–344.

Im, K. S., M. H. Pesaran, and Y. Shin. 2003. Testing for unit roots in heterogeneous panels. Journal of Econometrics 115: 53–74.

Martinez-Zarzoso, I., and A. Bengochea-Morancho. 2004. Pooled mean group estimation of an environmental kuznets curve for CO 2 . Economics Letters 82: 121–126.

Pesaran, M. H., Y. Shin, and R. P. Smith. 1997. Estimating long-run relationships in dynamic heterogeneous panels. DAE Working Papers Amalgamated Series 9721.

———. 1999. Pooled mean group estimation of dynamic heterogeneous panels. Journal of the American Statistical Association 94: 621–634.

Pesaran, M. H., and R. P. Smith. 1995. Estimating long-run relationships from dynamic heterogeneous panels. Journal of Econometrics 68: 79–113.

Phillips, P. C. B., and H. R. Moon. 2000. Nonstationary panel data analysis: An overview of some recent developments. Econometric Reviews 19: 263–286.

R的可视化能力非常强大，可以做出各种炫目的图表，但因为是命令行操作，入门门槛比较高；Excel的图表对于一般的统计分析可视化够用，掌握规律以后，也可以灵活组合创建商业级的图表，图形交互界面，上手快，但把图表做精美也不容易。

---