Implementing Marketing Analytics in R: Part 2
Feb 15, 2020 • 16 Minute Read
Introduction
Marketing and customer-related decisions are a top priority for every business. With the help of statistical modeling and analytics, it's possible to support decision makers and help them make strategic decisions based on data, not just instincts. The integration of statistical modeling with marketing strategy can also be referred to as marketing analytics.
In this series of two guides, we explore techniques for implementing marketing analytics in R.
In the previous guide, Part 1, we covered:
-
Customer Churn Prediction
-
RFM Analysis
In this guide, Part 2, we will cover:
-
Clustering
-
Sales Forecasting
-
Other Areas
Let's start by loading the required libraries.
library(plyr)
library(readr)
library(dplyr)
library(caret)
library(ggplot2)
library(repr)
library(caret)
Segmentation
Segmentation is the marketing strategy of dividing or segmenting customers to build and deploy effective marketing plans. The basis of segmenting can be demographic, geographic, or any other metric. The statistical technique of dividing the customers into homogeneous segments is called clustering. Let's start by loading the data.
Data
For this clustering exercise, we'll use a fictitious dataset of retail store customers containing 150 observations and 5 variables, as described below:
-
mean_HH_income: The mean household income of the customer.
-
ABV: Average bill value per purchase by a customer.
-
No_bills_month: Average number of bills in a month from the customer.
-
satisfaction_score: The customer's satisfaction score. The higher the satisfaction score, the more satisfied the customer is with the store.
-
family_size: Number of members in the customer's family.
-
Gender: Whether the customer is Female ('1') or male ('0').
Let's load the data and look at the structure.
df_seg = read_csv("segdata.csv")
glimpse(df_seg)
Output:
Observations: 150
Variables: 6
$ mean_HH_income <int> 88006, 88005, 88004, 88003, 88002, 88001, 88000, 64...
$ ABV <int> 355, 354, 353, 352, 351, 350, 349, 348, 347, 346, 3...
$ No_bills_month <int> 41, 21, 21, 53, 21, 44, 49, 22, 10, 58, 25, 52, 41,...
$ satisfaction_score <int> 1, 3, 4, 4, 2, 1, 3, 4, 4, 3, 5, 1, 5, 3, 3, 5, 4, ...
$ family_size <int> 4, 5, 3, 7, 6, 9, 6, 5, 8, 10, 3, 3, 3, 8, 7, 3, 3,...
$ Gender <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
The next step is to look at the summary of the numerical variables.
summary(df_seg)
Output:
mean_HH_income ABV No_bills_month satisfaction_score
Min. :45000 Min. : 26.00 Min. : 5.00 Min. :1.000
1st Qu.:52349 1st Qu.: 62.25 1st Qu.:16.00 1st Qu.:2.000
Median :64277 Median :286.50 Median :26.50 Median :3.000
Mean :59873 Mean :215.39 Mean :30.20 Mean :2.873
3rd Qu.:64314 3rd Qu.:323.75 3rd Qu.:44.75 3rd Qu.:4.000
Max. :88006 Max. :355.00 Max. :60.00 Max. :5.000
family_size Gender
Min. : 1.000 Min. :0.0000
1st Qu.: 3.000 1st Qu.:0.0000
Median : 6.000 Median :0.0000
Mean : 5.513 Mean :0.3333
3rd Qu.: 8.000 3rd Qu.:1.0000
Max. :10.000 Max. :1.0000
From the output above, we can infer that the variables have different scales that will adversely affect the calculation of distances for clustering. To prevent this, we'll normalize the data using the code below. We can normalize the variables in a data frame by using the preProcess function in the caret package.
The first line of code below pre-processes the data, while the second line performs the normalization. If we look at the summary of the segNorm variable, we'll see that all the variables now have a mean of zero.
p1_object = preProcess(df_seg)
segNorm = predict(p1_object, df_seg)
summary(segNorm)
Output:
mean_HH_income ABV No_bills_month satisfaction_score
Min. :-1.6354 Min. :-1.5124 Min. :-1.4713 Min. :-1.31892
1st Qu.:-0.8274 1st Qu.:-1.2229 1st Qu.:-0.8291 1st Qu.:-0.61487
Median : 0.4842 Median : 0.5679 Median :-0.2160 Median : 0.08918
Mean : 0.0000 Mean : 0.0000 Mean : 0.0000 Mean : 0.00000
3rd Qu.: 0.4883 3rd Qu.: 0.8654 3rd Qu.: 0.8495 3rd Qu.: 0.79323
Max. : 3.0935 Max. : 1.1149 Max. : 1.7399 Max. : 1.49728
family_size Gender
Min. :-1.5167 Min. :-0.7047
1st Qu.:-0.8446 1st Qu.:-0.7047
Median : 0.1635 Median :-0.7047
Mean : 0.0000 Mean : 0.0000
3rd Qu.: 0.8356 3rd Qu.: 1.4095
Max. : 1.5077 Max. : 1.4095
After normalizing the data, the next step is to run the clustering algorithm. We'll use the popular kmeans algorithm because it is powerful and interpretable. The number of clusters we will create is five, as given by the argument centers=5 in the first line of code. The second line displays the number of observations in each cluster.
kmeans_Clust = kmeans(segNorm, centers=5, iter.max=1000)
table(kmeans_Clust$cluster)
Output:
1 2 3 4 5
31 44 13 21 41
From the output above, we can infer that the smallest cluster, 3, has 13 observations, while the largest cluster, 2, has 44 observations. If we want to understand the clusters, we can do that with the command below.
The first line of code below adds a new variable in the original un-normalized data, df-seg. The second line uses the lapply() function to create and print the summary of each cluster.
df_seg$kmeansseg = kmeans_Clust$cluster
lapply(split(df_seg, df_seg$kmeansseg), colMeans)
Output:
$`1`
mean_HH_income ABV No_bills_month satisfaction_score
51832.129032 46.645161 22.709677 2.258065
family_size Gender kmeansseg
4.354839 1.000000 1.000000
$`2`
mean_HH_income ABV No_bills_month satisfaction_score
67507.954545 311.727273 36.818182 1.909091
family_size Gender kmeansseg
5.204545 0.000000 2.000000
$`3`
mean_HH_income ABV No_bills_month satisfaction_score
51437.230769 129.076923 28.692308 2.153846
family_size Gender kmeansseg
7.384615 0.000000 3.000000
$`4`
mean_HH_income ABV No_bills_month satisfaction_score
5.182395e+04 1.424762e+02 4.519048e+01 3.619048e+00
family_size Gender kmeansseg
6.142857e+00 9.047619e-01 4.000000e+00
$`5`
mean_HH_income ABV No_bills_month satisfaction_score
64556.170732 304.292683 21.560976 4.219512
family_size Gender kmeansseg
5.804878 0.000000 5.000000
Sales Forecasting
Forecasting needs no introduction, and it is an important part of marketing teams, often used for forecasting sales. Let's start by loading the data.
Data
For this forecasting exercise, we'll use fictitious monthly sales data containing 81 observations, and the task will be to create a forecasting model to predict sales over the next 12 months.
library(forecast)
library(fpp2)
library(TTR)
library(dplyr)
dat <- read_csv("salesforecasting.csv")
glimpse(dat)
Output:
Observations: 81
Variables: 3
$ Date <chr> "01-Mar-13", "01-Apr-13", "01-May-13", "01-Jun-13", "01-Jul-13",...
$ Sales <int> 192, 198, 197, 213, 202, 186, 215, 213, 213, 220, 210, 219, 233,...
$ Class <chr> "Train", "Train", "Train", "Train", "Train", "Train", "Train", "...
The next step is to create training and test sets for model validation. The first couple of lines of code below perform this task by using the Class variable, which already has labels Train and Test. The third line of code prints the number of rows in the training and the test data. Note that the test data set has 12 rows because we will be building the model to forecast for 12 months ahead.
dat_train = subset(dat, Class == 'Train')
dat_test = subset(dat, Class == 'Test')
nrow(dat_train); nrow(dat_test)
Output:
1] 69
[1] 12
Preparing the Time Series Object
To run forecasting models in 'R', we need to convert the data into a time series object, which is done in the first line of code below. The start and end argument specifies the time of the first and the last observations, respectively. The argument frequency specifies the number of observations per unit of time.
We will also create a utility function for calculating Mean Absolute Percentage Error (or MAPE), which will be used to evaluate the performance of the forecasting models. The lower the MAPE value, the better the forecasting model. This is done in the second to fourth lines of code.
dat_ts <- ts(dat_train[, 2], start = c(2013, 3), end = c(2018, 11), frequency = 12)
#lines 2 to 4
mape <- function(actual,pred){
mape <- mean(abs((actual - pred)/actual))*100
return (mape)
}
With the data and the mape function prepared, we move to the forecasting technique to be used in this case.
Holt's Trend Method
We'll use Holt's Trend forecasting method, a traditional statistical algorithm, which considers the trend component while generating forecasts. This method involves two smoothing equations, one for the level and one for the trend component.
The first line of code below creates the Holt's winter model and stores it in an object 'holt_model'. The second line prints the summary and the forecasts for the next 12 months.
holt_model <- holt(dat_ts, h = 12)
summary(holt_model)
Output:
Forecast method: Holt's method
Model Information:
Holt's method
Call:
holt(y = dat_ts, h = 12)
Smoothing parameters:
alpha = 0.8342
beta = 1e-04
Initial states:
l = 193.7099
b = 3.4287
sigma: 16.2391
AIC AICc BIC
682.6969 683.6492 693.8674
Error measures:
ME RMSE MAE MPE MAPE MASE ACF1
Training set -0.8870433 15.76137 11.8633 -0.4614957 4.084311 0.2366019 0.004106845
Forecasts:
Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
Dec 2018 382.4585 361.6473 403.2698 350.6305 414.2866
Jan 2019 385.8799 358.7765 412.9833 344.4288 427.3309
Feb 2019 389.3012 357.1116 421.4908 340.0715 438.5309
Mar 2019 392.7225 356.1461 429.2989 336.7838 448.6612
Apr 2019 396.1438 355.6521 436.6355 334.2171 458.0706
May 2019 399.5651 355.5036 443.6267 332.1789 466.9514
Jun 2019 402.9865 355.6225 450.3504 330.5496 475.4233
Jul 2019 406.4078 355.9563 456.8593 329.2489 483.5667
Aug 2019 409.8291 356.4676 463.1906 328.2197 491.4385
Sep 2019 413.2504 357.1288 469.3721 327.4198 499.0811
Oct 2019 416.6717 357.9187 475.4248 326.8168 506.5267
Nov 2019 420.0931 358.8209 481.3652 326.3854 513.8008
The output above shows that the MAPE for the training data is 4.1 percent. Let's now evaluate the model performance on the test data, which is done in the lines of code below. The MAPE error on the test data comes out to be 2.7 percent, which is an improvement over the training data result.
df_holt = as.data.frame(holt_model)
dat_test$holt = df_holt$`Point Forecast`
mape(dat_test$Sales, dat_test$holt)
Output:
1] 2.735082
Overall Forecasting for Next 12 Months
The accuracy number we got is respectable, so we can create the forecast for the next 12 months using the line of codes below.
df_forecast <- ts(dat[, 2], start = c(2013, 3), end = c(2019, 11), frequency = 12)
arima_model_ahead <- auto.arima(df_forecast)
forecast_12 = forecast::forecast(arima_model_ahead, h=12)
forecast_12
Output:
Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
Dec 2019 438.1473 419.0833 457.2113 408.9915 467.3031
Jan 2020 441.9584 417.3683 466.5485 404.3510 479.5658
Feb 2020 446.0598 416.9752 475.1444 401.5787 490.5409
Mar 2020 442.3232 409.3512 475.2952 391.8968 492.7496
Apr 2020 443.8119 407.3648 480.2591 388.0708 499.5530
May 2020 445.3007 405.6821 484.9193 384.7092 505.8921
Jun 2020 449.9827 407.4283 492.5370 384.9014 515.0640
Jul 2020 459.3094 414.0091 504.6097 390.0286 528.5902
Aug 2020 463.9914 416.1024 511.8804 390.7515 537.2313
Sep 2020 462.8675 412.5227 513.2122 385.8718 539.8631
Oct 2020 462.3241 409.6379 515.0103 381.7476 542.9007
Nov 2020 463.5226 408.5947 518.4505 379.5176 547.5275
Other Areas
Marketing analytics is a humongous area, and it is not possible to cover all the concepts and techniques in one guide, or for that matter even one book. We have covered the most widely used marketing analytics techniques above. Some other popular ones are mentioned below.
-
Customer Life Time Value (or CLV): In this technique, you can decide which customers are most valuable to the business by modeling their CLV.
-
Natural Language Processing (or NLP): NLP is a domain in itself and is increasingly being used for sentiment analysis, topic modeling, and text classification.
-
Dimensionality Reduction: This technique is often used for the survey of CRM data, which can be extensive.
-
Affinity Analysis: In this technique, you can decide which items have a higher likelihood of being bought together.
-
Marketing Campaign Optimization: In this technique, you can model the effectiveness of a marketing campaign and optimize your marketing budget.
Conclusion
In this series of guides, you learned about the most popular use cases of marketing analytics. You learned the concepts along with the implementation in R.
To learn more about data science with R, please refer to the following guides:
