TV-ad Attribution, Gaussian Processes

Problem description:

This work was done in Mister Spex GmbH, and slides of a presentation I gave at PyData meetup are available now here.

There is a website, in this case an e-shop, and we have information about user sessions on the website. We also have information about TV-ads shown to the public requested by the company. The question is, which of those sessions on the website are there because of the TV-ads.

There are some obstacles to answer the above question:

• Users don’t tell the website why they came to the website, or how they found about the service.
• Users use different channels to get to the website. They might have seen an ad on TV and then have used google to find the exact address to the website.

As a result, we change the question to the following: what is the likelihood of each session to be because of a TV-ad shown recently. Hence, the above input/output is given/desired:

• Input:
  • One entry for each session which are date-time labeled.
  • A list of date-time labeled events.
• Output:
  • What is the likelihood of each session to be influenced by a nearby event.

Prepare the data

We reduce the granularity of the data and aggregate it by minute. As a result, we have session count per minute, and the time of (TV-ad) events by minute. This results in a time series data having some points on it labeled by our events.

One way to see if the events have had an influence on our traffic, is to forget about those events, detect anomalies on our time series data, and see if they’ve happened around those events. I’ve used some of those methods in other projects, and they don’t make me happy on this data. They either detect wrong signals as anomalies, or miss some not so obvious signals.

Before going further, let’s look at the data. Fig. 1 illustrates a week worth of data. The X axis represents time (1 for each minute), and the Y axis shows the normalized session count. The data is normalized as X / median(X) simply to anonymize the actual session counts, since the actual numbers are irrelevant to this article.

Fig. 1: Normalized session count for a week

Some of those sharp peaks are right after a TV-ad event. Field knowledge says most of the traffic caused by a TV-ad comes up to 8 minutes after the event. I take precaution, and assume some error in reported time of the event. Now let’s remove any reported session from 2 minutes before the event, to 20 minutes after the event. Fig. 2 and 3 show one day of the data, before and after removing data around reported events.

Fig. 2: Normalized session count for a day

Fig. 3: Normalized session count for a day, after removing data around reported events

Method

Intuitively speaking, I would like to remove the parts of the data close to the given events, fit a model to the remaining data, and then see how close the observed data is to the model’s predictions. If the observed values are significantly higher than the predicted values, most probably some of those sessions are because of the event.

Given an observed and predicted values, deciding whether or not the observed value is significantly larger, can be tricky, and in some cases very arbitrary. Using a method which gives you the probability of an observed value at a given input would make the process much easier. One way of achieving such a goal is to have a method which, intuitively speaking, tells you what it thinks the predicted value is, and also how confident it is. A Gaussian Process (GP) is such a model, in a way that it gives a normal distribution at a given input point as the output. This predicted output has a mean and a variance, which we use to see how probable a given observation is, under the assumption that the model correctly describes the data.

Given a predicted normal distribution, the further the observed value is from the expected mean, the less expected it is. In other words, the larger the blue area shown in Fig. 4, the more probable it’s an outlier.

Fig. 4: Normal distribution and three different observed values - credit: thefemalecelebrity.com

The area under the curve marked by the blue areas shown in Fig. 4 is between 0 and 1. It is 0 if the expected value is the same as our observed value, and is close to 1 if it is the two observed and expected mean are far from each other. Hereafter score(x) represents this area under the curve. The piece of code bellow shows how this is calculated in python:

import math
import scipy

def phi(x):
    return 0.5 + 0.5 * scipy.special.erf(x / math.sqrt(2))

def score(x):
    return 1 - abs(phi(x) - phi(-x))

x_score = score((x - expected_mean) / expected_std)

Now for a moment, assume we know a given data point (number of sessions) is influenced by a TV-ad. Still not all those sessions are because of the event. Some of them would have been there even with no TV-ad event. We take the predicted mean as the background session count, and the rest as TV-ad influenced sessions. Therefore, if x is observed value, and x’ is predicted value, (x - x’) / x is the portion of sessions attributed to a TV-ad event.

At the end, we take score(x) * (x - x’) / x as the likelihood of each session being attributed to a TV-ad.

Implementation and Results

I used Python3, matplotlib to generate the plots, and GPy to train the Gaussian Process model.

There are at least two ways of fitting and using GPs in our setting. One is to train a really nice model using the historical data, and use the saved trained model to get predictions for future data. By looking at the data, I would sum these three kernels and optimize on the data:

• A periodic kernel to handle periodicity
• A Gaussian (RBF) kernel to handle the non-periodic part of the data
• A white noise kernel to handle fluctuations seen in the data

Because there were no big enough server available, training a model on the historical data at once was not feasible. Instead, for each day, I take a small part of the data which includes 10 hours before and after start and end of the day, and train a model on that part of the data. This makes training of the model on a small machine feasible, and it also doesn’t require the periodic kernel, which makes the optimization process even faster. Therefore I only use these two kernels:

• A Gaussian (RBF) kernel to handle the non-periodic part of the data
• A white noise kernel to handle fluctuations seen in the data

To summarize, here’s what I do:

• For a given date, take the 44 hours of data corresponding to the given date (10 hours before and 10 hours after the beginning and the end of the day respectively)
• Remove the part of the data around given any TV-ad events in that range.
• Fit a GP to the data.
• Calculate the abovementioned scores for each data point which corresponds to one minute.

Using GPy, I train the model as:

kernel = GPy.kern.RBF(input_dim=1) + GPy.kern.White(input_dim=1)
m = GPy.models.GPRegression(Xtr.reshape(-1,1),ytr.reshape(-1,1),kernel)
m.optimize()

In the above code, Xtr and ytr are the input/output vectors respectively. The Xtr is a sequence of integers representing the minute distance from where the data started, and ytr is the normalized number of sessions for that minute. This results in a model depicted in Fig. 5. The blue gradient represents the variance of the expected output, showing the expected value should be somewhere in the blue region. Please note that we’re only interested in the regions inside the data, and not the region outside the data range.

Fig. 5: The trained Gaussian Process model on the 44h of data.

Table 1 shows calculated scores and likelihoods of each session in each minute to be influenced by a TV-ad. A negative portion and likelihood simply mean the expected value is larger that the observed, and can safely be ignored. The Is Significant column simply flags if the score is over 0.9. I would personally attribute only sessions of those minutes to a TV-ad. Also, as expected, not all TV-ad events result in significant rise in the traffic, and some clearly result in more traffic than others.

Table 1: A piece of output of the method