Mixed Media Model
DW Richardson
2023-02-20
Establishing Environment
This R code chunk loads several R packages: datarium, PerformanceAnalytics, mctest, forecast, lmtest, car, and ggfortify. These packages are used for various tasks such as data exploration, statistical modeling, and visualization.
The datarium package provides tools for working with data sets in R, including functions for data import, cleaning, and manipulation.
PerformanceAnalytics provides a range of tools for performance and risk analysis, including functions for calculating and visualizing performance metrics, such as Sharpe ratios and maximum drawdowns.
The mctest package provides functions for performing multiple comparison tests in R, which can be used to compare the means of multiple groups or variables.
The forecast package provides tools for time series forecasting, including functions for modeling, fitting, and predicting time series data.
The lmtest package provides tools for hypothesis testing and diagnostics for linear regression models in R.
The car package provides tools for regression modeling and diagnostics, including functions for calculating influential observations and performing model comparisons.
The ggfortify package provides extensions to the ggplot2 package for creating visualizations of statistical models, including time series models and regression models.
library(datarium)
library(PerformanceAnalytics)## Loading required package: xts## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric##
## Attaching package: 'PerformanceAnalytics'## The following object is masked from 'package:graphics':
##
## legendlibrary(mctest)
library(forecast)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(lmtest)
library(car)## Loading required package: carDatalibrary(ggfortify)## Loading required package: ggplot2## Registered S3 methods overwritten by 'ggfortify':
## method from
## autoplot.Arima forecast
## autoplot.acf forecast
## autoplot.ar forecast
## autoplot.bats forecast
## autoplot.decomposed.ts forecast
## autoplot.ets forecast
## autoplot.forecast forecast
## autoplot.stl forecast
## autoplot.ts forecast
## fitted.ar forecast
## fortify.ts forecast
## residuals.ar forecastLoading Data
In the next three lines, the marketing
dataset is loaded from the datarium
package and assigned to the object
sampledf. The
str() function is used to display the
structure of the sampledf object, which
shows the variable names, data types, and number of observations. The
View() function is used to open a
spreadsheet-like view of the data for further examination.
data(marketing)
sampledf <- marketing
#View it
str(sampledf)## 'data.frame': 200 obs. of 4 variables:
## $ youtube : num 276.1 53.4 20.6 181.8 217 ...
## $ facebook : num 45.4 47.2 55.1 49.6 13 ...
## $ newspaper: num 83 54.1 83.2 70.2 70.1 ...
## $ sales : num 26.5 12.5 11.2 22.2 15.5 ...View(sampledf)Correlation Exploration
The chart.Correlation() function from
the PerformanceAnalytics package is used
to create a scatter plot matrix that displays the pairwise correlations
between the variables in a data frame. The function takes several
arguments, including the data frame to use, whether to show histograms
of the variables along the diagonal of the plot, and the symbol to use
for the points in the scatter plots.
In this case, the function is applied to the
sampledf data frame, which is a subset of
the marketing data included in the
datarium package. The
histogram argument is set to
TRUE, which displays histograms of each
variable in the diagonal plots, and the
pch argument is set to 19, which sets the
plotting symbol for the points in the scatter plots to a solid circle.
The resulting plot shows the pairwise correlations between the variables
in sampledf, as well as the distribution
of each variable.
chart.Correlation(sampledf, histogram = TRUE, pch=19)
Taking Stock of Adstock
This code chunk generates adstocked versions of the sample data frame’s Facebook, Youtube, and Newspaper columns. The adstocked versions of the columns represent the memory of previous weeks’ ad spending on each platform. First, it sets the adstock rate, which is the degree to which the effect of the ad spending decays over time. It sets the memory, which is the number of past weeks of ad spending that are included in the adstocked column. Then, it calculates the adstocked version of each column using a filter that applies the adstock formula to the column data. The resulting adstocked column data is then plotted against the original column data in a line plot, where the adstocked data is represented by the blue line. This allows for visual comparison of the impact of ad spending in each platform over time.
#set adstock fb rate
set_rate_fb <- 0.1
set_memory <- 2
get_adstock_fb <- rep(set_rate_fb, set_memory+1) ^ c(0:set_memory)
#set adstock youtube rate
set_rate_yt <- 0.15
set_memory <- 2
get_adstock_youtube <- rep(set_rate_yt, set_memory+1) ^ c(0:set_memory)
#set adstock news rate
set_rate_news <- 0.25
set_memory <- 2
get_adstock_news <- rep(set_rate_news, set_memory+1) ^ c(0:set_memory)
#adstocked fb
ads_fb <- stats::filter(c(rep(0, set_memory), sampledf$facebook), get_adstock_fb, method="convolution")
ads_fb <- ads_fb[!is.na(ads_fb)]
#plot
plot(seq(1,length(sampledf$facebook)), sampledf$facebook, type="h",
main = "Adstocked Facebook",
xlab="Time (Weeks)", ylab="Facebook",
ylim=c(0, max(c(sampledf$facebook, ads_fb))),
frame.plot=FALSE)
lines(ads_fb, col="blue")
#adstocked youtube
ads_youtube <- stats::filter(c(rep(0, set_memory), sampledf$youtube), get_adstock_youtube, method="convolution")
ads_youtube <- ads_youtube[!is.na(ads_youtube)]
#plot
plot(seq(1,length(sampledf$youtube)), sampledf$youtube, type="h",
main = "Adstocked Youtube",
xlab="Time (Weeks)", ylab="Youtube",
ylim=c(0, max(c(sampledf$youtube, ads_youtube))),
frame.plot=FALSE)
lines(ads_youtube, col="blue")
#adstocked newpaper
ads_news <- stats::filter(c(rep(0, set_memory), sampledf$newspaper), get_adstock_news, method="convolution")
ads_news <- ads_news[!is.na(ads_news)]
#plot
plot(seq(1,length(sampledf$newspaper)), sampledf$newspaper, type="h",
main = "Adstocked Newspaper",
xlab="Time (Weeks)", ylab="Newspaper",
ylim=c(0, max(c(sampledf$newspaper, ads_news))),
frame.plot=FALSE)
lines(ads_news, col="blue")
Initial Model & Diagnostics
The code performs multiple diagnostic tests and analyses on a multiple linear regression model.
The first line of code,
mmm_1 <- lm(sampledf$sales ~ ads_youtube + ads_fb + ads_news),
specifies a multiple linear regression model with sales as the dependent
variable and adstocked youtube, facebook, and newspaper as independent
variables. The lm() function is used to
fit the model, and the result is stored in
mmm_1.
The next line of code, summary(mmm_1),
provides a summary of the regression model’s coefficients, standard
errors, t-statistics, p-values, and other relevant information.
The imcdiag() function from the
mctest package is used to check for
multicollinearity. The method argument is
set to “VIF” to calculate variance inflation factors (VIFs) for each
independent variable. Alternatively, the
jtools package could be used to perform
the same analysis with the summ()
function.
The code then proceeds to check for heteroscedasticity. First, a 2x2
layout of residual plots is created with
par(mfrow=c(2,2)). The
plot() function is then used to create a
residuals vs fitted plot and a Scale-Location plot. These plots are used
to visually inspect the relationship between the fitted values and
residuals to detect any heteroscedasticity.
Finally, two formal tests for heteroscedasticity are conducted. The
bptest() function from the
lmtest package performs a Breusch-Pagan
test, while the ncvTest() function from
the car package performs a non-constant
variance test. These tests provide objective statistical evidence to
confirm the presence of heteroscedasticity in the model’s errors.
#We are specifying sales as the dependent variable in the lm() function
mmm_1 <- lm(sampledf$sales ~ ads_youtube + ads_fb + ads_news)
summary(mmm_1)##
## Call:
## lm(formula = sampledf$sales ~ ads_youtube + ads_fb + ads_news)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.0023 -1.2174 0.3924 1.4328 4.0541
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.947558 0.460560 4.229 3.6e-05 ***
## ads_youtube 0.045187 0.001479 30.557 < 2e-16 ***
## ads_fb 0.188177 0.009220 20.410 < 2e-16 ***
## ads_news -0.005953 0.006041 -0.985 0.326
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.168 on 196 degrees of freedom
## Multiple R-squared: 0.8819, Adjusted R-squared: 0.8801
## F-statistic: 487.9 on 3 and 196 DF, p-value: < 2.2e-16#Check for multicollinearity using VIFs
imcdiag(mmm_1, method = "VIF")##
## Call:
## imcdiag(mod = mmm_1, method = "VIF")
##
##
## VIF Multicollinearity Diagnostics
##
## VIF detection
## ads_youtube 1.0032 0
## ads_fb 1.1504 0
## ads_news 1.1497 0
##
## NOTE: VIF Method Failed to detect multicollinearity
##
##
## 0 --> COLLINEARITY is not detected by the test
##
## ===================================#check for heteroscedasticity
#first, plot the model out and review the siduals vs fitted plot and the Sclae-Location plot
par(mfrow=c(2,2)) # put all 4 charts into 1 page
plot(mmm_1)
#Confirm with an objective test for heteroscedasticity using Breusch Pagan test and NCV test
lmtest::bptest(mmm_1)##
## studentized Breusch-Pagan test
##
## data: mmm_1
## BP = 4.1583, df = 3, p-value = 0.2449car::ncvTest(mmm_1)## Non-constant Variance Score Test
## Variance formula: ~ fitted.values
## Chisquare = 2.063199, Df = 1, p = 0.15089Taking Time Into Account
In this R code chunk, the time series analysis and forecasting are performed. Firstly, a time series is created using the sales data, which is a weekly data set, with the frequency of 52. It is decomposed to get the individual components for trends, seasonality, etc. The time series linear regression model, tslm(), is used for regression. In the model, the trend and season values are automatically generated based on the ts() object. The model mmm_2 is created and a summary of the model is printed out.
After building the model, the focus is on forecasting. A data frame is created that contains new values for advertising spending. 40% of the newspaper spending is assigned to YouTube, 60% to Facebook, and the newspaper spending is set to zero. Then all new values are put into a dataframe. Finally, the model is used to predict sales for the next period based on these new budget allocation values.
#frequency is 52 to denote weekly as there are about 52 weeks in a year.
#ts() needs a minimum of 2 periods (52 x 2 = 104 weeks),
#our data has observations from 200 weeks so this should be sufficient
ts_sales <- ts(sampledf$sales, start = 1, frequency = 52)
#check class. should state "ts"
class(ts_sales)## [1] "ts"#decompose to get the individual components for trends, seasonality, etc
ts_sales_comp <- decompose(ts_sales)
#plot out
plot(ts_sales_comp)
#we use tslm() for our regression.
#this is just a timeseries wrapper for lm() but allows trend and seasons on the fly from the data
#https://www.rdocumentation.org/packages/forecast/versions/8.16/topics/tslm
#just specify "trend" and "season" and tslm() will automatically generate values based on the ts() object you have specified.
#fit the model
mmm_2 <- tslm(ts_sales ~ trend + season + ads_youtube + ads_fb + ads_news)
summary(mmm_2)##
## Call:
## tslm(formula = ts_sales ~ trend + season + ads_youtube + ads_fb +
## ads_news)
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.8701 -1.0978 0.1012 1.1498 5.5102
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.550710 1.269482 3.585 0.000461 ***
## trend -0.001346 0.002698 -0.499 0.618585
## season2 -2.439879 1.500236 -1.626 0.106066
## season3 -3.837931 1.511244 -2.540 0.012160 *
## season4 -0.019777 1.505982 -0.013 0.989540
## season5 -2.589598 1.533250 -1.689 0.093391 .
## season6 -2.975709 1.500929 -1.983 0.049317 *
## season7 -1.828409 1.498110 -1.220 0.224279
## season8 -1.737503 1.499661 -1.159 0.248538
## season9 -1.735511 1.542630 -1.125 0.262446
## season10 -2.379873 1.520560 -1.565 0.119747
## season11 -3.564521 1.504577 -2.369 0.019158 *
## season12 -2.221964 1.503207 -1.478 0.141552
## season13 -2.126054 1.502541 -1.415 0.159235
## season14 -1.757011 1.537121 -1.143 0.254913
## season15 -1.008243 1.507831 -0.669 0.504776
## season16 -0.971702 1.504299 -0.646 0.519340
## season17 -1.675260 1.510448 -1.109 0.269230
## season18 -1.487987 1.499527 -0.992 0.322714
## season19 -3.228179 1.516448 -2.129 0.034975 *
## season20 -0.920816 1.496122 -0.615 0.539217
## season21 -2.694057 1.496636 -1.800 0.073942 .
## season22 -2.724598 1.534971 -1.775 0.078008 .
## season23 -4.770156 1.513552 -3.152 0.001976 **
## season24 -2.457553 1.519227 -1.618 0.107930
## season25 -0.131696 1.520747 -0.087 0.931110
## season26 -2.774732 1.523675 -1.821 0.070671 .
## season27 -5.345331 1.525779 -3.503 0.000612 ***
## season28 -1.525467 1.516993 -1.006 0.316301
## season29 -3.807714 1.501743 -2.536 0.012295 *
## season30 -1.881424 1.501888 -1.253 0.212343
## season31 -2.258311 1.511899 -1.494 0.137444
## season32 -2.855228 1.502809 -1.900 0.059442 .
## season33 -2.696570 1.502977 -1.794 0.074887 .
## season34 -2.413012 1.510219 -1.598 0.112282
## season35 -2.534971 1.521820 -1.666 0.097937 .
## season36 -2.044690 1.504609 -1.359 0.176287
## season37 -0.511780 1.509610 -0.339 0.735092
## season38 -1.460185 1.502723 -0.972 0.332833
## season39 -1.605796 1.502239 -1.069 0.286887
## season40 0.044458 1.530591 0.029 0.976868
## season41 -0.503530 1.515917 -0.332 0.740250
## season42 -1.542739 1.511224 -1.021 0.309035
## season43 -2.440774 1.510239 -1.616 0.108250
## season44 -2.802112 1.506946 -1.859 0.065002 .
## season45 -4.075788 1.628560 -2.503 0.013443 *
## season46 -1.791812 1.624069 -1.103 0.271743
## season47 -1.542609 1.628470 -0.947 0.345085
## season48 -1.246402 1.616600 -0.771 0.441969
## season49 -3.781911 1.643492 -2.301 0.022819 *
## season50 -1.074279 1.637860 -0.656 0.512933
## season51 -4.164205 1.645888 -2.530 0.012480 *
## season52 -3.191509 1.648129 -1.936 0.054772 .
## ads_youtube 0.044778 0.001669 26.825 < 2e-16 ***
## ads_fb 0.185387 0.010483 17.685 < 2e-16 ***
## ads_news -0.008845 0.006874 -1.287 0.200297
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.11 on 144 degrees of freedom
## Multiple R-squared: 0.9178, Adjusted R-squared: 0.8865
## F-statistic: 29.25 on 55 and 144 DF, p-value: < 2.2e-16#we want to forecast using our model
#i.e. if we were to spend x1 in youtube, x2 in facebook and 0 in newspaper for the next period what would this look like?
#we first need to create a dataframe containing new figures
#we want to get newspaper spend
news_spend <- as.data.frame(ads_news)
names(news_spend)[1] <- "ads_news"
#and give 40% to Youtube.
#this is added to the current youtube spend budget (assuming we are keeping the youtube budget
#for the next period is the same as it has been for the previous period)
yt_spend <- as.data.frame(ads_youtube)
names(yt_spend)[1] <- "ads_youtube"
yt_spend$ads_youtube <- yt_spend$ads_youtube + (news_spend$ads_news*0.4)
yt_spend## ads_youtube
## 1 309.3360
## 2 124.7700
## 3 75.6147
## 4 223.8465
## 5 281.8254
## 6 91.8375
## 7 97.4796
## 8 165.4629
## 9 36.0855
## 10 255.1974
## 11 129.7062
## 12 280.3926
## 13 101.8287
## 14 138.5649
## 15 288.0336
## 16 304.9785
## 17 184.4907
## 18 397.2108
## 19 154.4226
## 20 209.8518
## 21 318.9354
## 22 346.4301
## 23 92.6988
## 24 301.9788
## 25 129.6264
## 26 345.2001
## 27 229.4211
## 28 334.0293
## 29 359.7543
## 30 159.0057
## 31 397.2246
## 32 215.0442
## 33 165.1983
## 34 344.1663
## 35 169.7604
## 36 378.2142
## 37 378.8319
## 38 168.3219
## 39 94.8543
## 40 304.3179
## 41 305.2647
## 42 278.3340
## 43 396.1035
## 44 319.9560
## 45 99.2952
## 46 241.3323
## 47 162.0507
## 48 322.8627
## 49 345.4869
## 50 151.8603
## 51 280.4574
## 52 165.2343
## 53 303.6246
## 54 293.8188
## 55 369.8148
## 56 323.3652
## 57 79.2039
## 58 184.8603
## 59 299.0031
## 60 303.9474
## 61 120.3366
## 62 357.9819
## 63 355.9485
## 64 182.3181
## 65 197.9661
## 66 113.3229
## 67 55.7907
## 68 179.8800
## 69 317.3745
## 70 321.3351
## 71 306.5238
## 72 194.1276
## 73 71.5287
## 74 181.3596
## 75 290.7186
## 76 107.6088
## 77 62.8608
## 78 161.9883
## 79 35.7495
## 80 156.0675
## 81 126.4638
## 82 325.7250
## 83 156.2838
## 84 124.2036
## 85 292.0161
## 86 308.7768
## 87 148.6785
## 88 186.0174
## 89 171.2421
## 90 186.0189
## 91 196.1421
## 92 79.9566
## 93 302.5731
## 94 383.8122
## 95 195.5979
## 96 250.9353
## 97 278.9208
## 98 274.7121
## 99 413.6502
## 100 248.2143
## 101 329.9859
## 102 455.1474
## 103 419.4708
## 104 298.1058
## 105 332.5614
## 106 242.9223
## 107 82.7484
## 108 133.1733
## 109 48.6300
## 110 317.6868
## 111 345.8337
## 112 355.6608
## 113 266.0736
## 114 295.7919
## 115 154.2279
## 116 139.6122
## 117 202.3164
## 118 130.5207
## 119 208.9104
## 120 68.6208
## 121 203.6739
## 122 78.9228
## 123 290.9211
## 124 197.8836
## 125 341.1780
## 126 170.9817
## 127 60.8745
## 128 111.2634
## 129 283.1646
## 130 134.6874
## 131 26.9601
## 132 342.9522
## 133 64.2639
## 134 295.6224
## 135 121.0338
## 136 83.8416
## 137 47.8623
## 138 364.3791
## 139 118.8402
## 140 242.0769
## 141 129.5340
## 142 288.5313
## 143 329.0988
## 144 193.7679
## 145 164.1585
## 146 198.5202
## 147 322.3944
## 148 361.4241
## 149 107.1477
## 150 79.6914
## 151 366.5010
## 152 225.3669
## 153 280.2489
## 154 265.6560
## 155 271.0392
## 156 48.3561
## 157 143.6976
## 158 214.6677
## 159 69.6663
## 160 186.9516
## 161 251.2659
## 162 165.8319
## 163 265.2885
## 164 240.5289
## 165 179.4048
## 166 348.4845
## 167 87.5604
## 168 272.1615
## 169 326.8113
## 170 396.0816
## 171 128.3178
## 172 239.2281
## 173 68.8800
## 174 219.6555
## 175 306.0552
## 176 398.8788
## 177 369.0798
## 178 277.0143
## 179 385.5918
## 180 265.4694
## 181 232.0059
## 182 309.5352
## 183 128.7912
## 184 399.9855
## 185 381.7524
## 186 314.6112
## 187 227.1726
## 188 272.4810
## 189 386.1225
## 190 91.3017
## 191 64.1910
## 192 102.4929
## 193 51.3585
## 194 210.9945
## 195 214.3884
## 196 84.7416
## 197 129.6819
## 198 234.8454
## 199 407.5104
## 200 346.6590#and give the remainder 60% of newspaper spend to facebook.
#this is added to the current fb spend budget (assuming we are keeping the fb budget
#for the next period is the same as it has been for the previous period)
fb_spend <- as.data.frame(ads_fb)
names(fb_spend)[1] <- "ads_fb"
fb_spend$ads_fb <- fb_spend$ads_fb + (news_spend$ads_news*0.6)
fb_spend## ads_fb
## 1 95.1840
## 2 96.6240
## 3 121.3776
## 4 112.1631
## 5 74.1633
## 6 127.6161
## 7 78.4056
## 8 43.9998
## 9 9.1311
## 10 19.5732
## 11 28.5822
## 12 37.7172
## 13 94.3266
## 14 30.8460
## 15 78.1947
## 16 107.9712
## 17 143.7108
## 18 115.5609
## 19 58.1412
## 20 51.1722
## 21 79.0635
## 22 37.1223
## 23 62.3694
## 24 51.0987
## 25 37.4628
## 26 24.4278
## 27 49.1367
## 28 43.2315
## 29 56.0526
## 30 57.1809
## 31 75.6837
## 32 61.8720
## 33 34.7196
## 34 31.7418
## 35 10.8300
## 36 12.7935
## 37 58.5318
## 38 98.7717
## 39 72.2166
## 40 80.4513
## 41 61.6959
## 42 78.2004
## 43 47.1996
## 44 34.8783
## 45 68.1894
## 46 61.8468
## 47 48.2109
## 48 72.4215
## 49 64.9233
## 50 52.7445
## 51 39.0951
## 52 22.5084
## 53 81.9462
## 54 110.1132
## 55 64.4004
## 56 111.9939
## 57 81.3171
## 58 49.1088
## 59 94.1562
## 60 55.8114
## 61 25.3137
## 62 95.4885
## 63 54.2130
## 64 51.3159
## 65 78.6465
## 66 22.8792
## 67 34.1961
## 68 28.2441
## 69 44.8902
## 70 78.1770
## 71 75.5730
## 72 52.3728
## 73 63.0267
## 74 38.4081
## 75 46.5345
## 76 123.5949
## 77 39.0447
## 78 52.8894
## 79 49.5747
## 80 32.1330
## 81 53.9598
## 82 39.8379
## 83 56.2179
## 84 89.0277
## 85 89.3901
## 86 82.7640
## 87 60.5910
## 88 103.5813
## 89 100.7460
## 90 113.9712
## 91 31.1730
## 92 30.7086
## 93 89.2773
## 94 111.9990
## 95 45.0990
## 96 83.3415
## 97 22.4205
## 98 45.2817
## 99 94.4115
## 100 98.6220
## 101 57.0936
## 102 128.2539
## 103 52.3386
## 104 43.5681
## 105 51.3462
## 106 104.2419
## 107 51.4221
## 108 26.9418
## 109 24.5925
## 110 41.9436
## 111 55.8948
## 112 74.0283
## 113 31.5849
## 114 36.2040
## 115 85.6908
## 116 92.4987
## 117 51.3921
## 118 20.7315
## 119 105.3876
## 120 54.6156
## 121 75.3648
## 122 75.0555
## 123 28.1886
## 124 56.0724
## 125 99.2988
## 126 51.0132
## 127 92.9166
## 128 21.7071
## 129 65.5038
## 130 52.3020
## 131 63.7140
## 132 42.8415
## 133 43.1067
## 134 78.2838
## 135 106.1109
## 136 81.3915
## 137 64.0812
## 138 84.9645
## 139 60.9405
## 140 63.7353
## 141 36.4953
## 142 101.8773
## 143 85.7685
## 144 46.2408
## 145 54.7479
## 146 19.1544
## 147 18.8001
## 148 93.5658
## 149 71.2611
## 150 55.3515
## 151 51.1431
## 152 54.7086
## 153 49.7898
## 154 82.4283
## 155 44.6286
## 156 24.4389
## 157 91.6587
## 158 33.7617
## 159 84.1485
## 160 60.6651
## 161 54.7368
## 162 87.9318
## 163 54.9207
## 164 58.9161
## 165 28.6452
## 166 68.6466
## 167 76.7634
## 168 32.4648
## 169 75.3312
## 170 31.4634
## 171 32.4672
## 172 64.3272
## 173 48.3672
## 174 25.5918
## 175 17.6742
## 176 92.2032
## 177 64.8783
## 178 44.4498
## 179 28.3719
## 180 30.8916
## 181 14.5581
## 182 28.9260
## 183 34.2087
## 184 110.6238
## 185 66.6489
## 186 79.9350
## 187 32.2176
## 188 54.0072
## 189 27.2862
## 190 34.8654
## 191 59.4933
## 192 24.4542
## 193 30.8022
## 194 59.5716
## 195 54.1992
## 196 20.3940
## 197 15.3372
## 198 18.4794
## 199 100.7553
## 200 33.9396#leaving nothing for newspapers (we are swtiching it off for the next period)
final_news_spend <- as.data.frame(news_spend*0)
names(final_news_spend)[1] <- "ads_news"
final_news_spend## ads_news
## 1 0
## 2 0
## 3 0
## 4 0
## 5 0
## 6 0
## 7 0
## 8 0
## 9 0
## 10 0
## 11 0
## 12 0
## 13 0
## 14 0
## 15 0
## 16 0
## 17 0
## 18 0
## 19 0
## 20 0
## 21 0
## 22 0
## 23 0
## 24 0
## 25 0
## 26 0
## 27 0
## 28 0
## 29 0
## 30 0
## 31 0
## 32 0
## 33 0
## 34 0
## 35 0
## 36 0
## 37 0
## 38 0
## 39 0
## 40 0
## 41 0
## 42 0
## 43 0
## 44 0
## 45 0
## 46 0
## 47 0
## 48 0
## 49 0
## 50 0
## 51 0
## 52 0
## 53 0
## 54 0
## 55 0
## 56 0
## 57 0
## 58 0
## 59 0
## 60 0
## 61 0
## 62 0
## 63 0
## 64 0
## 65 0
## 66 0
## 67 0
## 68 0
## 69 0
## 70 0
## 71 0
## 72 0
## 73 0
## 74 0
## 75 0
## 76 0
## 77 0
## 78 0
## 79 0
## 80 0
## 81 0
## 82 0
## 83 0
## 84 0
## 85 0
## 86 0
## 87 0
## 88 0
## 89 0
## 90 0
## 91 0
## 92 0
## 93 0
## 94 0
## 95 0
## 96 0
## 97 0
## 98 0
## 99 0
## 100 0
## 101 0
## 102 0
## 103 0
## 104 0
## 105 0
## 106 0
## 107 0
## 108 0
## 109 0
## 110 0
## 111 0
## 112 0
## 113 0
## 114 0
## 115 0
## 116 0
## 117 0
## 118 0
## 119 0
## 120 0
## 121 0
## 122 0
## 123 0
## 124 0
## 125 0
## 126 0
## 127 0
## 128 0
## 129 0
## 130 0
## 131 0
## 132 0
## 133 0
## 134 0
## 135 0
## 136 0
## 137 0
## 138 0
## 139 0
## 140 0
## 141 0
## 142 0
## 143 0
## 144 0
## 145 0
## 146 0
## 147 0
## 148 0
## 149 0
## 150 0
## 151 0
## 152 0
## 153 0
## 154 0
## 155 0
## 156 0
## 157 0
## 158 0
## 159 0
## 160 0
## 161 0
## 162 0
## 163 0
## 164 0
## 165 0
## 166 0
## 167 0
## 168 0
## 169 0
## 170 0
## 171 0
## 172 0
## 173 0
## 174 0
## 175 0
## 176 0
## 177 0
## 178 0
## 179 0
## 180 0
## 181 0
## 182 0
## 183 0
## 184 0
## 185 0
## 186 0
## 187 0
## 188 0
## 189 0
## 190 0
## 191 0
## 192 0
## 193 0
## 194 0
## 195 0
## 196 0
## 197 0
## 198 0
## 199 0
## 200 0#now put these new values all into a dataframe.
#We'll use the model to predict sales for the next period based on these new budget allocation values
new_spends <- cbind(yt_spend, fb_spend,final_news_spend)
new_spends## ads_youtube ads_fb ads_news
## 1 309.3360 95.1840 0
## 2 124.7700 96.6240 0
## 3 75.6147 121.3776 0
## 4 223.8465 112.1631 0
## 5 281.8254 74.1633 0
## 6 91.8375 127.6161 0
## 7 97.4796 78.4056 0
## 8 165.4629 43.9998 0
## 9 36.0855 9.1311 0
## 10 255.1974 19.5732 0
## 11 129.7062 28.5822 0
## 12 280.3926 37.7172 0
## 13 101.8287 94.3266 0
## 14 138.5649 30.8460 0
## 15 288.0336 78.1947 0
## 16 304.9785 107.9712 0
## 17 184.4907 143.7108 0
## 18 397.2108 115.5609 0
## 19 154.4226 58.1412 0
## 20 209.8518 51.1722 0
## 21 318.9354 79.0635 0
## 22 346.4301 37.1223 0
## 23 92.6988 62.3694 0
## 24 301.9788 51.0987 0
## 25 129.6264 37.4628 0
## 26 345.2001 24.4278 0
## 27 229.4211 49.1367 0
## 28 334.0293 43.2315 0
## 29 359.7543 56.0526 0
## 30 159.0057 57.1809 0
## 31 397.2246 75.6837 0
## 32 215.0442 61.8720 0
## 33 165.1983 34.7196 0
## 34 344.1663 31.7418 0
## 35 169.7604 10.8300 0
## 36 378.2142 12.7935 0
## 37 378.8319 58.5318 0
## 38 168.3219 98.7717 0
## 39 94.8543 72.2166 0
## 40 304.3179 80.4513 0
## 41 305.2647 61.6959 0
## 42 278.3340 78.2004 0
## 43 396.1035 47.1996 0
## 44 319.9560 34.8783 0
## 45 99.2952 68.1894 0
## 46 241.3323 61.8468 0
## 47 162.0507 48.2109 0
## 48 322.8627 72.4215 0
## 49 345.4869 64.9233 0
## 50 151.8603 52.7445 0
## 51 280.4574 39.0951 0
## 52 165.2343 22.5084 0
## 53 303.6246 81.9462 0
## 54 293.8188 110.1132 0
## 55 369.8148 64.4004 0
## 56 323.3652 111.9939 0
## 57 79.2039 81.3171 0
## 58 184.8603 49.1088 0
## 59 299.0031 94.1562 0
## 60 303.9474 55.8114 0
## 61 120.3366 25.3137 0
## 62 357.9819 95.4885 0
## 63 355.9485 54.2130 0
## 64 182.3181 51.3159 0
## 65 197.9661 78.6465 0
## 66 113.3229 22.8792 0
## 67 55.7907 34.1961 0
## 68 179.8800 28.2441 0
## 69 317.3745 44.8902 0
## 70 321.3351 78.1770 0
## 71 306.5238 75.5730 0
## 72 194.1276 52.3728 0
## 73 71.5287 63.0267 0
## 74 181.3596 38.4081 0
## 75 290.7186 46.5345 0
## 76 107.6088 123.5949 0
## 77 62.8608 39.0447 0
## 78 161.9883 52.8894 0
## 79 35.7495 49.5747 0
## 80 156.0675 32.1330 0
## 81 126.4638 53.9598 0
## 82 325.7250 39.8379 0
## 83 156.2838 56.2179 0
## 84 124.2036 89.0277 0
## 85 292.0161 89.3901 0
## 86 308.7768 82.7640 0
## 87 148.6785 60.5910 0
## 88 186.0174 103.5813 0
## 89 171.2421 100.7460 0
## 90 186.0189 113.9712 0
## 91 196.1421 31.1730 0
## 92 79.9566 30.7086 0
## 93 302.5731 89.2773 0
## 94 383.8122 111.9990 0
## 95 195.5979 45.0990 0
## 96 250.9353 83.3415 0
## 97 278.9208 22.4205 0
## 98 274.7121 45.2817 0
## 99 413.6502 94.4115 0
## 100 248.2143 98.6220 0
## 101 329.9859 57.0936 0
## 102 455.1474 128.2539 0
## 103 419.4708 52.3386 0
## 104 298.1058 43.5681 0
## 105 332.5614 51.3462 0
## 106 242.9223 104.2419 0
## 107 82.7484 51.4221 0
## 108 133.1733 26.9418 0
## 109 48.6300 24.5925 0
## 110 317.6868 41.9436 0
## 111 345.8337 55.8948 0
## 112 355.6608 74.0283 0
## 113 266.0736 31.5849 0
## 114 295.7919 36.2040 0
## 115 154.2279 85.6908 0
## 116 139.6122 92.4987 0
## 117 202.3164 51.3921 0
## 118 130.5207 20.7315 0
## 119 208.9104 105.3876 0
## 120 68.6208 54.6156 0
## 121 203.6739 75.3648 0
## 122 78.9228 75.0555 0
## 123 290.9211 28.1886 0
## 124 197.8836 56.0724 0
## 125 341.1780 99.2988 0
## 126 170.9817 51.0132 0
## 127 60.8745 92.9166 0
## 128 111.2634 21.7071 0
## 129 283.1646 65.5038 0
## 130 134.6874 52.3020 0
## 131 26.9601 63.7140 0
## 132 342.9522 42.8415 0
## 133 64.2639 43.1067 0
## 134 295.6224 78.2838 0
## 135 121.0338 106.1109 0
## 136 83.8416 81.3915 0
## 137 47.8623 64.0812 0
## 138 364.3791 84.9645 0
## 139 118.8402 60.9405 0
## 140 242.0769 63.7353 0
## 141 129.5340 36.4953 0
## 142 288.5313 101.8773 0
## 143 329.0988 85.7685 0
## 144 193.7679 46.2408 0
## 145 164.1585 54.7479 0
## 146 198.5202 19.1544 0
## 147 322.3944 18.8001 0
## 148 361.4241 93.5658 0
## 149 107.1477 71.2611 0
## 150 79.6914 55.3515 0
## 151 366.5010 51.1431 0
## 152 225.3669 54.7086 0
## 153 280.2489 49.7898 0
## 154 265.6560 82.4283 0
## 155 271.0392 44.6286 0
## 156 48.3561 24.4389 0
## 157 143.6976 91.6587 0
## 158 214.6677 33.7617 0
## 159 69.6663 84.1485 0
## 160 186.9516 60.6651 0
## 161 251.2659 54.7368 0
## 162 165.8319 87.9318 0
## 163 265.2885 54.9207 0
## 164 240.5289 58.9161 0
## 165 179.4048 28.6452 0
## 166 348.4845 68.6466 0
## 167 87.5604 76.7634 0
## 168 272.1615 32.4648 0
## 169 326.8113 75.3312 0
## 170 396.0816 31.4634 0
## 171 128.3178 32.4672 0
## 172 239.2281 64.3272 0
## 173 68.8800 48.3672 0
## 174 219.6555 25.5918 0
## 175 306.0552 17.6742 0
## 176 398.8788 92.2032 0
## 177 369.0798 64.8783 0
## 178 277.0143 44.4498 0
## 179 385.5918 28.3719 0
## 180 265.4694 30.8916 0
## 181 232.0059 14.5581 0
## 182 309.5352 28.9260 0
## 183 128.7912 34.2087 0
## 184 399.9855 110.6238 0
## 185 381.7524 66.6489 0
## 186 314.6112 79.9350 0
## 187 227.1726 32.2176 0
## 188 272.4810 54.0072 0
## 189 386.1225 27.2862 0
## 190 91.3017 34.8654 0
## 191 64.1910 59.4933 0
## 192 102.4929 24.4542 0
## 193 51.3585 30.8022 0
## 194 210.9945 59.5716 0
## 195 214.3884 54.1992 0
## 196 84.7416 20.3940 0
## 197 129.6819 15.3372 0
## 198 234.8454 18.4794 0
## 199 407.5104 100.7553 0
## 200 346.6590 33.9396 0Forecasting & Adjusting
This code chunk starts by resetting the number of plots per page to 1
with the command par(mfrow=c(1,1)) and
setting a random seed value.
The first step is to generate a forecast with no changes in the
budget allocations using the forecast()
function and assigning the output to the object
forecast_unchanged. The forecast is then
plotted using ggplot2::autoplot(), where
the time series is represented by a black line, the predicted values are
in red, and the shaded area represents the confidence interval.
Next, a new forecast is generated using the same model, but with new
data for the ad spends specified in
new_spends. This is accomplished by
including the argument newdata=new_spends
in the forecast() function. The forecast
is stored in the object
forecast_new_spends and is plotted in blue
using ggplot2::autoplot(), with the same
format as the first plot.
The two forecasts are then overlaid on the same plot using
forecast::autolayer(), which takes as
input the forecast_new_spends object and
specifies the color of the overlay. The resulting plot shows both
forecasts with the same time series and confidence intervals, but with
different predicted values in red and blue.
Finally, the fitted values from the second forecast can be accessed
using forecast_new_spends$fitted.
par(mfrow=c(1,1)) # reset to 1 chart per page
set.seed(9999)
#what performance looks like with no change
forecast_unchanged <- forecast(mmm_2, h=200)
ggplot2::autoplot(forecast_unchanged, ts.colour = 'black', size= 0.7, predict.size = 0.7, predict.colour = 'red', conf.int = TRUE, conf.int.fill = 'red', main = "Forecasted", predict.linetype='dashed') 
#forecast with budget changes
forecast_new_spends <- forecast(mmm_2, newdata=new_spends)
ggplot2::autoplot(forecast_new_spends, ts.colour = 'black', size= 0.7, predict.size = 0.7, predict.colour = 'blue', conf.int = TRUE, conf.int.fill = 'blue', main = "Forecasted")
#overlaying them together using autolayer()
forecast_unchanged <- forecast(mmm_2, h=200)
ggplot2::autoplot(forecast_unchanged, ts.colour = 'black', size= 0.7, predict.size = 0.7, predict.colour = 'red', conf.int = TRUE, conf.int.fill = 'red', main = "Forecasted", predict.linetype='dashed') + forecast::autolayer(forecast_new_spends, col = 'blue')
#Get fitted values
#Finally, you can access fitted values from the model by quering forecast_new_spends$fitted
forecast_new_spends$fitted## Time Series:
## Start = c(1, 1)
## End = c(4, 44)
## Frequency = 52
## 1 2 3 4 5 6 7 8
## 24.588115 15.275364 12.538218 22.320221 15.513745 13.895686 14.012805 14.698709
## 9 10 11 12 13 14 15 16
## 5.222714 13.541318 7.182988 19.976709 13.125628 10.736301 22.319081 26.491980
## 17 18 19 20 21 22 23 24
## 16.213091 27.741890 12.470354 17.988623 20.944636 17.891381 5.724418 18.421034
## 25 26 27 28 29 30 31 32
## 12.506749 17.469268 15.468322 21.443906 22.334114 12.408564 24.989036 14.149712
## 33 34 35 36 37 38 39 40
## 8.630951 21.696425 10.029641 20.048003 30.547049 21.025241 12.772188 25.857938
## 41 42 43 44 45 46 47 48
## 22.178255 21.864505 26.311774 17.635804 9.220648 17.687279 11.475277 26.255609
## 49 50 51 52 53 54 55 56
## 18.810911 11.636555 12.354779 10.451529 26.202482 24.211837 23.664639 28.496526
## 57 58 59 60 61 62 63 64
## 11.002516 13.795223 26.128716 23.416866 8.495647 25.779695 19.884029 16.832566
## 65 66 67 68 69 70 71 72
## 20.366005 10.588397 11.589339 14.988835 23.049515 26.785333 21.314693 14.802221
## 73 74 75 76 77 78 79 80
## 11.737037 10.657902 17.659616 14.052656 7.104476 14.596040 7.531763 11.501844
## 81 82 83 84 85 86 87 88
## 11.554080 17.257277 12.460942 16.026440 23.980477 18.545055 14.069995 18.112354
## 89 90 91 92 93 94 95 96
## 15.328150 20.224016 12.956954 7.363324 22.789529 26.075124 13.692829 18.347198
## 97 98 99 100 101 102 103 104
## 13.710652 18.888036 29.450094 22.619873 15.317170 28.252094 20.497926 17.750602
## 105 106 107 108 109 110 111 112
## 26.979590 21.996567 6.391892 9.661761 3.044299 21.217480 18.572526 26.074639
## 113 114 115 116 117 118 119 120
## 18.490185 19.862829 17.442850 15.306908 14.076247 8.266823 18.282338 9.490861
## 121 122 123 124 125 126 127 128
## 16.411452 9.857985 14.261926 19.498223 22.413089 11.281260 9.342249 7.041108
## 129 130 131 132 133 134 135 136
## 27.634025 9.963678 8.790097 18.220087 9.215163 22.239176 14.450755 15.765750
## 137 138 139 140 141 142 143 144
## 13.121990 23.660935 12.485488 23.224969 13.996541 21.636026 23.850344 13.602809
## 145 146 147 148 149 150 151 152
## 13.094650 11.795111 17.617723 27.290305 14.388701 12.024684 21.594629 13.524518
## 153 154 155 156 157 158 159 160
## 17.591919 23.111351 17.307295 6.117869 19.069813 11.596232 10.205252 16.081492
## 161 162 163 164 165 166 167 168
## 16.159440 15.291611 17.965953 21.489786 14.471454 15.936158 11.770132 15.323817
## 169 170 171 172 173 174 175 176
## 20.232120 22.488479 11.169241 17.308325 9.725942 14.094792 15.193026 30.750934
## 177 178 179 180 181 182 183 184
## 24.945238 15.129457 16.473717 15.205202 14.914750 15.691015 4.929817 27.794164
## 185 186 187 188 189 190 191 192
## 21.696643 26.094983 12.659267 19.338099 22.306582 8.217586 13.694875 9.894674
## 193 194 195 196 197 198 199 200
## 6.248260 21.354717 20.860514 9.335929 10.457567 15.145261 27.457674 19.046693