MAPAx example for R

By | November 9, 2017

I have been asked a few time to provide an example how to use MAPAx, from the MAPA package for R, so I prepared this blog post. I admit the documentation could be better, so I put together this example from a retailing case – the original setting MAPAx was developed for (see paper here). This provides a step by step use of MAPAx.

First we load the package and some data.


y <- cbind(c(1030.829,893.551,1084.09,1278.436,936.708,915.322,885.713,2364.399,774.88,977.506,831.616,813.656,1569.956,967.925,806.146,1063.117,1123.787,906.686,996.498,1088.464,977.414,1128.328,896.594,1007.172,1046.379,1514.648,1626.115,2959.558,838.506,949.377,1433.307,805.048,1218.907,872.43,1730.103,865.734,1845.713,919.291,1003.363,1102.969,847.38,1965.26,809.673,953.193,1066.089,991.352,1115.694,1003.333,1090.48,930.749,1006.184,1239.068,873.707,728.583,881.316,1302.468,997.442,3481.118,841.042,997.601,1830.194,909.693,2358.406,2573.673,777.08,773.781,945.424,968.646,1074.589,1046.22,1155.559,990.627,931.943,786.285,2297.025,628.166,889.238,937.631,1113.925,870.384,1018.375,799.458,1542.328,1879.587,750.307,1087.948,1247.803,1052.352,883.899,793.126,913.736,1082.142,968.823,2099.176,841.224,964.227,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA))
x <- cbind(c(0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0),
 c(0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0),
 c(0,0,0,1,1,0,0,1,1,1,1,1,1,1,0,0,0,0,0,0,0,1,1,1,1,1,0,0,0,0,1,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,1,1,1,0,1,0,0,0,1,1,1,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,1,0,0,1,1,1,0,0,1,1,0,0))
y <- y[!is.na(y)]
y <- ts(y,frequency=12)

This toy set includes sales of a heavily promoted item, with three different promotion. Let’s visualise what we have:


plot(y)
cmp <- brewer.pal(3,"Set1")
for (i in 1:3){
  points(time(y)[which(x[,i]==1)],y[which(x[,i]==1)],col=cmp[i],pch=i,cex=1.5)
}
legend("topright",c("Promo1","Promo2","Promo3"),col=cmp,pch=1:3)

Now let’s fit a MAPA and MAPAx and produce forecasts. I am modelling the sales in logs to capture the multiplicative promotional effects. For other types of applications this is not needed.

library(MAPA)
mapafit.x <- mapaest(log(y),type="es",display=1,outplot=1,xreg=x)
mapafit <- mapaest(log(y),type="es",outplot=1)

frc.x <- mapafor(log(y),mapafit.x,ifh=13,fh=13,xreg=x,conf.lvl=c(0.8,0.9,0.95),comb="w.mean")
frc <- mapafor(log(y),mapafit,ifh=13,fh=13,conf.lvl=c(0.8,0.9,0.95),comb="w.mean")

Let’s plot the results


par(mfrow=c(1,2))
plot(1:96, y, type="l",xlim=c(1,109),xaxs="i",main="MAPA",xlab="Period",ylab="Sales")
for (i in 1:96){
  lines(i:(12+i),exp(frc$infor[,i]),col="red")
}
cmp <- brewer.pal(9,"Reds")[4:2]
for (i in 1:3){
  polygon(c(97:109,109:97),exp(c(frc$PI[i,],rev(frc$PI[7-i,]))),col=cmp[i],border=NA)
}
lines(97:109,exp(frc$outfor),col="red")

plot(1:96, y, type="l",xlim=c(1,109),xaxs="i",main="MAPAx",xlab="Period",ylab="Sales")
for (i in 1:96){
  lines(i:(12+i),exp(frc.x$infor[,i]),col="red")
}
cmp <- brewer.pal(9,"Reds")[4:2]
for (i in 1:3){
  polygon(c(97:109,109:97),exp(c(frc.x$PI[i,],rev(frc.x$PI[7-i,]))),col=cmp[i],border=NA)
}
lines(97:109,exp(frc.x$outfor),col="red")

As you can see, MAPA provides a flat forecast. This is “correct” in the sense that apart from promotions this time series contains like information that could be captured by exponential smoothing. On the other hand, MAPAx is provided with the promotional information and makes use of it.

10 thoughts on “MAPAx example for R

  1. AnscombesGimlet

    Really appreciate the how-to, but can you fix the formatting? Also, I think you may have forgot the library(MAPA) line (should be obvious, but you never know if a new R user is trying the code). Thanks again for the example!

    Reply
  2. forecaster

    Hi Nikos,

    I have this issue when I select “MAPA” package. i know this is due to “smooth” package but raised it because I use only MAPA and smooth is required for “MAPA”. I have raised it as an issue in github for forecast package. I have cross posted it below. Appreciate all your help to the forecasting community.

    I get following error when I use forecast function. The problem arises due to conflict with “smooth” package. I need to select smooth package because “MAPA” package requires it.

    Error in UseMethod(“forecast”) :
    no applicable method for ‘forecast’ applied to an object of class “c(‘ARIMA’, ‘Arima’)”

    Reply
  3. Pedro

    Hi Nikos,

    When I run
    mapafit_x = mapaest(y,type=”es”,display=1,outplot=1,xreg=as.matrix(x))

    I get this error
    Aggregation level: 80/365 (21.92%)Error in if (results[[besti]]$ICs[ic] <= results[[i]]$ICs[ic]) { : missing value where TRUE/FALSE needed

    Do you know what may be the reason?

    Reply
    1. Nikos Post author

      I suspect the time series is very short, but I cannot be sure.

      Reply
  4. Danny

    Hello Nikos,

    Thanks a lot for making MAPA available!
    Just like in your example, I’m using 3 regressors and the function mapafit.x is working fine.
    However, when I run frc.x I get the following:

    “Error in mapacalc(inobs, mapafit, fh = ifh.c, comb, outplot = 0, hybrid, :
    Number of xreg input variables does not match mapafit specification.”

    Do you have an hypothesis for why this is happening? Many thanks in advance

    Reply
    1. Nikos Post author

      Not sure. The easiest way to get the xreg part working is to imput the complete matrix in both the estimation and forecast functions. It should be able to find what sample to use on its own from the target variable. It does not use dates to do that, but it counts observations, so the X’s do not have to be time series.

      Reply
  5. Ioannis

    Hi Nikos

    MAPAx sounds like an interesting alternative that would be worth trying out.

    My experience with using ARIMAX and ETS so far in predicting 7-day forecasts for SKUs is that both models are very sensitive to trends just before the forecasting horizon: e.g. if the last 1-2 days see an increase in sales, both models assume an increase in every single day of the horizon leading to bad forecasting accuracy especially in the latest 3-7 days of the horizon where they extrapolate that initial trend too much .

    If you have such experience yourself, do you know how MAPAx deals with ‘last-minute’ trends that do not last?

    Thanks,
    Ioannis

    Reply
    1. Nikos Post author

      By construction MAPA and MAPAx would be resistant to such changes, as short term movements would be filtered out in the temporally aggregate views of the time series. Therefore, even if on the original sampling frequency the ETS model might have poor performance, at an aggregate level it would be more robust and eventually the final MAPA forecast that would combine information from both would be resistant to these effects.

      Reply

Leave a Reply

Your email address will not be published. Required fields are marked *