#install and load package arules#install.packages("arules")
library(arules)#install and load arulesViz#install.packages("arulesViz")
library(arulesViz)#install and load tidyverse#install.packages("tidyverse")
library(tidyverse)#install and load readxml#install.packages("readxml")
library(readxl)#install and load knitr#install.packages("knitr")
library(knitr)#load ggplot2 as it comes in tidyverse
library(ggplot2)#install and load lubridate#install.packages("lubridate")
library(lubridate)#install and load plyr#install.packages("plyr")
library(plyr)
library(dplyr)技术支持
数据预处理
使用 read_excel(path to file) 将下载的文件读入 R。请在 read_excel(path-to-file-with-filename) 中填写包含文件名的完整路径。
#read excel into R dataframe
retail <- read_excel('D:/Documents/Online_Retail.xlsx')#complete.cases(data) will return a logical vector indicating which rows have no missing values. Then use the vector to get only rows that are complete using retail[,].
retail <- retail[complete.cases(retail),]#mutate function is from dplyr package. It is used to edit or add new columns to dataframe. Here Description column is being converted to factor column. as.factor converts column to factor column. %>% is an operator with which you may pipe values to another function or expression
retail %>% mutate(Description = as.factor(Description))技术支持
#Converts character data to date. Store InvoiceDate as date in new variable
retail$Date <- as.Date(retail$InvoiceDate)#Extract time from InvoiceDate and store in another variable
TransTime<- format(retail$InvoiceDate,"%H:%M:%S")#Convert and edit InvoiceNo into numeric
InvoiceNo <- as.numeric(as.character(retail$InvoiceNo))技术支持
您需要做的是,将 retail 数据框按 CustomerID、CustomerID 与 Date 分组,或者按 InvoiceNo 与 Date 分组。我们需要进行这样的分组,然后对其应用一个函数,并将输出存入另一个数据框。这可以通过 ddply 完成。
下面的代码将把同一 InvoiceNo 与 date 的所有产品合并,将该 InvoiceNo 与 date 下的所有商品合并为一行,各商品之间用 , 分隔。
library(plyr)#ddply(dataframe, variables_to_be_used_to_split_data_frame, function_to_be_applied)
transactionData <- ddply(retail,c("InvoiceNo","Date"),function(df1)paste(df1$Description,
collapse =","))#The R function paste() concatenates vectors to character and separated results using collapse=[any optional charcater string ]. Here ',' is used技术支持
write.csv(transactionData,"D:/Documents/market_basket_transactions.csv", quote =FALSE, row.names =FALSE)#transactionData: Data to be written#"D:/Documents/market_basket.csv": location of file with file name to be written to#quote: If TRUE it will surround character or factor column with double quotes. If FALSE nothing will be quoted#row.names: either a logical value indicating whether the row names of x are to be written along with x, or a character vector of row names to be written.技术支持
检查一下您的交易数据格式是否正确:
接下来,需要将该交易数据加载为 transaction 类的对象。这可以使用 arules 包中的 R 函数 read.transactions 完成。
tr <- read.transactions('D:/Documents/market_basket_transactions.csv', format ='basket', sep=',')#sep tell how items are separated. In this case you have separated using ','技术支持
transactions in sparse format with
22191 transactions (rows) and
7876 items (columns)
transactions as itemMatrix in sparse format with
22191 rows (elements/itemsets/transactions) and
7876 columns (items) and a density of 0.001930725
most frequent items:
WHITE HANGING HEART T-LIGHT HOLDER REGENCY CAKESTAND 3 TIER
18031709
JUMBO BAG RED RETROSPOT PARTY BUNTING
14601285
ASSORTED COLOUR BIRD ORNAMENT (Other)1250329938
element (itemset/transaction) length distribution:
sizes
12345678910111213141516173598159411419088617586966766635936245375165315515224641819202122232425262728293031323334441483419395315306272238253229213222215170159138142353637383940414243444546474849505113410911190113949387886563676360594964525354555657585960616263646566676840414943362939302728172525202724226970717273747576777879808182838485152019131616111512791415128911868788899091929394959697989910010110211148656116443652424103104105106107108109110111112113114116117118120121432263432131333121221231251261271311321331341401411421431451461471502132211211221121115415716817117717818020220422823624925028532040041932221111111111111
Min. 1st Qu. Median Mean 3rd Qu. Max.
1.003.0010.0015.2121.00419.00
includes extended item information - examples:
labels
11 HANGER
210 COLOUR SPACEBOY PEN
312 COLOURED PARTY BALLOONS技术支持
# Create an item frequency plot for the top 20 itemsif(!require("RColorBrewer")){# install color package of R
install.packages("RColorBrewer")#include library RColorBrewer
library(RColorBrewer)}
itemFrequencyPlot(tr,topN=20,type="absolute",col=brewer.pal(8,'Pastel2'), main="Absolute Item Frequency Plot")技术支持
set of 49122 rules
rule length distribution (lhs + rhs):sizes
234567891010521116854164241485561021937613121
Min. 1st Qu. Median Mean 3rd Qu. Max.
2.0005.0005.0005.4996.00010.000
summary of quality measures:
support confidence lift count
Min. :0.001036 Min. :0.8000 Min. :9.846 Min. :23.001st Qu.:0.0010821st Qu.:0.83331st Qu.:22.2371st Qu.:24.00
Median :0.001262 Median :0.8788 Median :28.760 Median :28.00
Mean :0.001417 Mean :0.8849 Mean :64.589 Mean :31.453rd Qu.:0.0015323rd Qu.:0.92593rd Qu.:69.2003rd Qu.:34.00
Max. :0.015997 Max. :1.0000 Max. :715.839 Max. :355.00
mining info:
data ntransactions support confidence
tr 221910.0010.8技术支持
# Here lhs=METAL because you want to find out the probability of that in how many customers buy METAL along with other items
inspect(head(metal.association.rules))技术支持
Apriori
Parameter specification:
confidence minval smax arem aval originalSupport maxtime support minlen maxlen target
0.80.11 none FALSETRUE50.001110 rules
ext
FALSE
Algorithmic control:
filter tree heap memopt load sort verbose
0.1TRUETRUEFALSETRUE2TRUE
Absolute minimum support count:22
set item appearances ...[1 item(s)] done [0.00s].
set transactions ...[30066 item(s),22191 transaction(s)] done [0.10s].
sorting and recoding items ...[2324 item(s)] done [0.02s].
creating transaction tree ... done [0.02s].
checking subsets of size 12 done [0.01s].
writing ...[1 rule(s)] done [0.00s].
creating S4 object ... done [0.01s].技术支持
# Here lhs=METAL because you want to find out the probability of that in how many customers buy METAL along with other items
inspect(head(metal.association.rules))技术支持
lhs rhs support confidence lift count
[1]{METAL}=>{DECORATION}0.0022531661443.8250技术支持
# Filter rules with confidence greater than 0.4 or 40%
subRules<-association.rules[quality(association.rules)$confidence>0.4]#Plot SubRules
plot(subRules)技术支持
'plotly_arules' is deprecated.
Use 'plot' instead.
See help("Deprecated")plot: Too many rules supplied. Only plotting the best 1000 rules using measure lift (change parameter max if needed)To reduce overplotting, jitter is added! Use jitter =0 to prevent jitter.技术支持
您是一名数据科学家(或正在成为数据科学家!),接到了一位零售店主的委托。客户向您提供了一段时间内店内多位顾客的所有交易数据,包含所购买的商品,并请您利用这些数据来助力业务增长。客户将根据您的发现,不仅调整、更新或新增库存商品,还可能据此调整线下门店或线上商店的布局。为得到能帮助客户的结果,您将使用市场篮子分析(MBA),它会在给定的交易数据上进行关联规则挖掘。
在本教程中,您将学习:
关联规则挖掘
当您希望在一个集合中寻找不同对象之间的关联、在交易数据库、关系型数据库或其他信息库中寻找频繁模式时,就会用到关联规则挖掘。关联规则挖掘的应用包括市场营销、零售中的篮子数据分析(又称市场篮子分析)、聚类与分类。它可以通过生成一组称为关联规则的规则,告诉您顾客经常一起购买哪些商品。简单来说,它给出的输出是“如果…则…”形式的规则。客户可以将这些规则用于多种营销策略:
请看下面的示例:
给定一组交易数据。可以看到交易编号 1 到 5。每笔交易显示了该交易中购买的商品。您可以看到,尿布与啤酒在三笔交易中同时出现。类似地,面包与牛奶在三笔交易中同时出现,使它们都成为频繁项集。关联规则以如下形式给出:
$A=>B[Support,Confidence]$
$=>$ 之前的部分称为if(前件,Antecedent),$=>$ 之后的部分称为then(后件,Consequent)。
其中 A 和 B 是交易数据中的商品集合。A 与 B 互不相交。
$Computer=>Anti-virus Software[Support=20\%,confidence=60\%]$
上述规则表示:
接下来一节,您将了解关联规则挖掘的基本概念:
关联规则挖掘的基本概念
项集(Itemset):由一个或多个商品组成的集合。K-项集表示包含 k 个商品的集合。
支持计数(Support Count):某个项集出现的频次
支持度(s):包含项集 X 的交易占所有交易的比例
$Support(X)=\frac{frequency(X)}{N}$
对于规则 A=>B,支持度为:
$Support(A=>B)=\frac{frequency(A,B)}{N}$
注:P(A∪B) 表示 A 与 B 同时发生的概率。P 表示概率。
请动手试试,计算 Milk=>Diaper 的支持度。
$Confidence(A=>B)=\frac{P(A\cap B)}{P(A)}=\frac{frequency(A,B)}{frequency(A)}$
即同时包含 A 和 B 的交易数,除以包含 A 的交易总数。
$Confidence(Bread=>Milk)=\frac{3}{4}=0.75=75\%$
现在请计算 Milk=>Diaper 的置信度。
注:支持度与置信度衡量规则的有趣程度。它们由最小支持度与最小置信度阈值设定。由客户设定的这些阈值有助于按自身需求比较规则强度。越接近阈值,规则对客户的价值通常越高。
频繁项集:支持度大于或等于最小支持度阈值(min_sup)的项集。在上例中 min_sup=3。该值由用户自行设定。
强规则:若规则 A=>B[Support, Confidence] 同时满足 min_sup 与 min_confidence,则为强规则。
提升度(Lift):在规则 A=>B 中,提升度刻画 A 与 B 的相关性。相关性反映项集 A 对项集 B 的影响程度。
$Lift(A=>B)=\frac{Support}{Supp(A)Supp(B)}$
例如,对于规则
{Bread}=>{Milk},提升度计算如下:$support(Bread)=\frac{4}{5}=0.8$
$support(Milk)=\frac{4}{5}=0.8$
$Lift(Bread=>Milk)=\frac{0.6}{0.8*0.8}=0.9$
若规则的提升度为 1,则 A 与 B 相互独立,无法从中推导规则。
若提升度 > 1,则 A 与 B 存在依赖关系,其强度由提升度数值体现。
若提升度 < 1,则 A 的出现会对 B 产生负向影响。
关联规则挖掘的目标
当您对给定的交易集合 T 应用关联规则挖掘时,目标是找到所有满足以下条件的规则:
APRIORI 算法
在本教程的这一部分,您将了解 R 库在进行市场篮子分析时所依赖的算法。这将帮助您更好地理解客户,并以更专注的方式开展分析。如果您已经了解 APRIORI 算法及其原理,可以直接跳到编码部分。
关联规则挖掘可视为两步法:
频繁项集生成:找出所有支持度 ≥ 预设最小支持计数的频繁项集
规则生成:从频繁项集中列出所有关联规则,计算各规则的支持度与置信度,剪除不满足最小支持度与最小置信度阈值的规则。
频繁项集生成是计算量最大的一步,因为它需要完整扫描数据库。
在上述步骤中,频繁项集生成的计算开销最大。
前面示例中只有 5 笔交易,但在现实世界中,零售的交易数据可能达到 GB 乃至 TB 级别。为在后续步骤前剪除无用项集,需要一个优化算法。这就是 APRIORI 算法。它指出:
该原理可用项集格(Itemset Lattice)来表示,这是对 APRIORI 原理的图形化展示。它由 k-项集的节点以及该 k-项集子集之间的关系构成。
如上图所示,底部是交易数据中的所有商品,然后向上逐层创建子集,直至空集。若有 d 个商品,格的规模为 $2^d$。这说明为每种组合计算支持度以生成频繁项集将有多难。下图展示了 APRIORI 如何显著减少需要生成的集合数量:
如果项集 {a,b} 不频繁,那么它的所有超集都无需考虑。
我们通过一个示例来理解。下面的示例中,您将看到为什么 APRIORI 是一种高效算法,并将逐步生成强关联规则。请拿起纸笔,跟着做!
如您所见,首先为 1-项集创建候选列表,其中包含交易数据中出现的所有单个商品。设想真实世界的零售交易数据,您就会明白候选生成有多昂贵。此处 APRIORI 发挥作用,帮助减少候选列表数量,最终生成有用规则。接下来的步骤将展示我们如何完成关联规则挖掘第一步——频繁项集生成。
下一步是列出所有频繁项集。您将取最后一个非空的频繁项集,本例中为 L2={I1, I2},{I2, I3}。然后对该频繁项集列表中的每个项集,生成其所有非空子集。请按下图所示继续:
可以看到,上面共有四条强规则。例如,${I2}=>I3$ 的置信度为 75%,表示购买了 I2 的人中有 75% 同时也购买了 I3。
至此,您已经学习了完整的 APRIORI 算法——它是数据挖掘中使用最广的算法之一。下面进入代码部分,呼——!
使用 R 实现 MBA/关联规则挖掘
本教程将使用 UCI 机器学习资料库中的一个数据集。该数据集名为 Online-Retail,您可以从这里下载。该数据集包含一家总部在英国、无实体店的在线零售商在 2010/12/01 至 2011/12/09 的交易数据。之所以使用这个数据集而非 R 自带数据集,是因为您更可能收到这种形式的零售数据,需要先进行数据预处理。
数据集描述
属性信息
加载库
首先,加载所需的库。下面的表格给出了这些库的简要说明(摘自这里),帮助您了解各库的作用:
arulesarulesViztidyversereadxlplyrggplot2knitrlubridate数据预处理
使用
read_excel(path to file)将下载的文件读入 R。请在read_excel(path-to-file-with-filename)中填写包含文件名的完整路径。现在,数据框(dataframe) retail 将包含 10 个属性,新增了
Date与Time两列。在应用 MBA/关联规则挖掘前,我们需要将数据框转换为交易数据,使得同一发票中一起购买的所有商品位于同一行。您可从
glimpse的输出中看到,每笔交易都是原子形态,即同一发票下的所有商品在关系型数据库中以原子记录表示。这种格式也称为singles 格式。您需要做的是,将
retail数据框按CustomerID、CustomerID与Date分组,或者按InvoiceNo与Date分组。我们需要进行这样的分组,然后对其应用一个函数,并将输出存入另一个数据框。这可以通过ddply完成。下面的代码将把同一
InvoiceNo与date的所有产品合并,将该InvoiceNo与date下的所有商品合并为一行,各商品之间用,分隔。接下来,由于 InvoiceNo 和 Date 在规则挖掘中并无用处,您可以将它们设为
NULL。这种交易数据格式称为basket 格式。接下来,您需将该交易数据存为一个 .csv(逗号分隔值)文件。使用
write.csv()即可。检查一下您的交易数据格式是否正确:
接下来,需要将该交易数据加载为 transaction 类的对象。这可以使用
arules包中的 R 函数read.transactions完成。下面的代码会将 basket 格式的交易数据文件
D:/Documents/market_basket_transactions.csv转换为 transaction 类对象。运行上述代码时,输出中可能出现大量
EOF within quoted string,无需担心。如果您已有数据框形式的交易数据,可用下行代码将其转换为 transaction 对象:
查看
tr这个 transaction 对象:summary(tr)是一个非常有用的命令,它为我们提供关于 transaction 对象的信息。来看一下上述输出的含义:共有 22191 笔交易(行)和 7876 个商品(列)。注意,7876 是数据集中涉及的产品描述数,22191 笔交易是这些商品的组合。
稀疏度中的密度(Density)表示稀疏矩阵中非零单元的百分比。也可理解为购买的商品总数,除以矩阵中可能的商品数量。您可用密度计算购买的商品数:22191x7876x0.001930725=337445
Summary 还可以告诉您最常见的商品。
元素(项集/交易)长度分布:展示了 1-项集、2-项集等对应的交易数量。第一行表示项个数,第二行表示交易数量。
例如,仅有 3598 笔只包含 1 个商品的交易,1594 笔包含 2 个商品的交易,而最长的一笔交易包含 419 个商品。
您可以生成
itemFrequencyPlot来创建条形图,基于itemMatrix(如>transactions或>itemsets、>rules中的商品)展示对象的频次分布,这正符合我们的场景。在
中,第一个参数是要绘制的交易对象
tr。topN用于绘制出现频次最高的前 N 个商品。type可为type="absolute"或type="relative"。若为 absolute,将绘制每个商品的绝对频次;若为 relative,将绘制这些商品相对其他商品出现的频率。该图显示,“WHITE HANGING HEART T-LIGHT HOLDER”和“REGENCY CAKESTAND 3 TIER”的销量最高。因此,为提高“SET OF 3 CAKE TINS PANTRY DESIGN”的销量,零售商可以将其摆放在“REGENCY CAKESTAND 3 TIER”附近。
您可以在这里探索
itemFrequencyPlot的其他选项。生成规则!
下一步是使用 APRIORI 算法挖掘规则。函数
apriori()来自arules包。apriori将tr作为要挖掘的交易对象。parameter用于设置最小支持度与最小置信度。其默认值为最小支持度 0.1、最小置信度 0.8、规则中最多 10 个商品(maxlen)。summary(association.rules)显示如下内容:参数设定:min_sup=0.001、min_confidence=0.8,规则中最多 10 个商品。
规则总数:共有 49122 条规则
规则长度分布:包含 5 个商品的规则最多:16424 条,2 个商品的规则最少:105 条
质量指标汇总:支持度、置信度与提升度的最小与最大值。
用于生成规则的信息:我们提供给算法的数据、支持度与置信度。
由于共有 49122 条规则,我们只打印前 10 条:
根据上述输出,您可以做出如下分析:
购买“WOBBLY CHICKEN”的顾客有 100% 也购买了“METAL”。
购买“BLACK TEA”的顾客有 100% 也购买了“SUGAR JARS”。
限制规则的数量与长度
如何限制生成规则的数量与长度?您可以在
apriori中设置parameters。通过调整这些参数来控制规则数量。若想得到更强的规则,可提高conf;若想得到更长的规则,可提高maxlen。移除冗余规则
您可以移除那些为更大规则子集的规则。使用下面的代码来移除这类规则:
which()返回向量中值为 TRUE 的元素位置。colSums()对数据框和数值数组计算行、列和。is.subset()判断一个向量的元素是否包含另一个向量的所有元素。查找与指定商品相关的规则
有时您想聚焦某个特定产品。若想找出影响商品 X 购买的因素,可以在
apriori命令中使用appearance选项。appearance允许我们设置规则的前件(IF 部分,LHS)和后件(THEN 部分,RHS)。例如,要找出顾客在购买“METAL”之前会买什么,请运行以下代码:
类似地,要回答“购买 METAL 的顾客还会购买……”这一问题,您需要将 METAL 放在 lhs:
可视化关联规则
由于根据数据会生成数百乃至数千条规则,您需要多种方式来展示发现。上文已介绍过的
ItemFrequencyPlot也是获取热销商品的好方法。这里将讨论以下可视化方式:
散点图
一种直接的关联规则可视化方式,是使用
arulesViz包的plot()绘制散点图。坐标轴使用支持度与置信度,此外默认用第三个指标提升度来为点着色(灰度)。上图显示提升度高的规则往往支持度较低。您可以为图形使用以下选项:
rulesObject:要绘制的规则对象measure:规则有趣度的度量。可以是支持度、置信度、提升度,或根据method的不同组合使用。shading:用于着色的度量(支持度、置信度、提升度)。默认是提升度。method:可视化方法(scatterplot、two-key plot、matrix3D)。two-key plot 在 x、y 轴上分别使用支持度和置信度,并用order(规则中的商品数量)来着色。
交互式散点图
一个很棒的交互式图可以用来展示您的规则,结合
arulesViz与plotly。您可以将鼠标悬停在每条规则上,查看所有质量指标(支持度、置信度与提升度)。基于图的可视化
基于图的技术使用顶点与边来可视化关联规则:顶点用商品名称标注,项集或规则作为第二组顶点表示。通过有向箭头将商品与项集/规则连接。由商品指向规则顶点的箭头表示 LHS 项,而从规则指向商品的箭头表示 RHS。顶点的大小与颜色通常代表兴趣度量。
图形绘制是直观的规则可视化方式,但当规则数量增多时容易拥挤。因此更适合用它来可视化较少数量的规则。
让我们从
subRules中选取置信度最高的 10 条规则。现在绘制一个交互式图:
注意:您可以在
plot中通过engine=htmlwidget参数将所有图形变为交互式。在
arulesViz中,关联规则集的图可以导出为 GraphML 格式或 Graphviz dot 文件,以便在 Gephi 等工具中探索。例如,导出提升度最高的 1000 条规则:单条规则表示
这种表示也称为平行坐标图。它有助于直观看到哪些产品与哪些商品组合会带来怎样的销售。
如上所述,RHS 为后件,即我们推测客户会购买的商品;LHS 的位置中,2 表示购物篮中最新加入的商品,1 表示先前已有的商品。
看最上方的箭头。它表明当我的购物车中已有“CHILDS GARDEN SPADE PINK”和“CHILDS GARDEN RAKE PINK”时,我也很可能同时购买“CHILDS GARDEN RAKE BLUE”。
结论
恭喜您!您已经学会了 APRIORI——数据挖掘中使用最频繁的算法之一。您了解了关联规则挖掘、其应用,以及它在零售中的应用——市场篮子分析。您也已经能够在 R 中实现市场篮子分析,并用精美的图形展示您的关联规则!继续加油学习!
参考文献:
https://datascienceplus.com/a-gentle-introduction-on-market-basket-analysis%E2%80%8A-%E2%80%8Aassociation-rules/
https://en.wikipedia.org/wiki/Sparse_matrix
https://cran.r-project.org/web/packages/arulesViz/vignettes/arulesViz.pdf
如果您想进一步学习 R,请参加 DataCamp 的 Importing and Managing Financial Data in R 课程。