Good Data Visualization Should Be Boring

So I’m going to make a statement that I’m sure some people are going to disagree with: good data visualization should be boring.

Well, at least kind of boring.

I’ve had a lot of conversations with a lot of people over the last few years or so about data visualization: why it’s important, what constitutes good and bad, and examples of its application in both problematic and very effective ways.

A salient point someone made to me once is that part of the problem with the practice of data visualization is that it isn’t viewed as a standalone discipline; it’s simply done, in high school math classes, university courses, or even in the workplace by professionals, and usually assumed that people will just pick it up without discussion around it and its proper application.

I think this is gradually starting to change, as with all the talk (or hype, depending on your point of view) around “Big Data”, analytics is becoming more mainstream, and data visualization is as well as a part of it. I also think dataviz is beginning to – gradually, very gradually – become viewed as a standalone discipline, with courses now being offered in it, and the “data visualization evangelism” of academics such as Edward Tufte and Alberto Cairo and work of practitioners like Stephen Few and Mike Bostock helping to raise awareness of what’s doing it wrong and what’s doing it right. This, along with others creating visualizations which go viral or delivering inspirational TED talks, are doing a lot for visualization as a practice.

The thing I found when I first started to get into dataviz is that even if you’re good with data that doesn’t necessarily mean you’re good at visualizing it. This is because, in addition to working with data, doing proper visualization involves questions of design and also the psychology of perception.

Less is More

I’m a minimalist, and therefore take what I call a functionalist perspective of data visualization. That is to say, the purpose of visualization is to most effectively represent that data so that it can be understood by the audience both most quickly and easily.

As such, I feel that good data visualization should be somewhat dull, or at least somewhat dry; in terms of depicting information and people perceiving it, it is usually the case that simpler is better. This is illustrated in principles like Tufte’s data-ink ratio.

So, look at the charts below. Which is more visually appealing to you? Which is simpler? Which one depicts the quantities such that you are able to interpret them the most quickly, accurately and with the most clarity?

If you’re like me, you’ll say the one on the right, which is a better visualization, even though it may not be as visually appealing to some. Most often you’re better served by a simpler, cleaner visualization (or perhaps several of them) than a lot of complexity and visual noise that doesn’t add to the reader’s understanding.

Never say always

That being said, as I mentioned, choices around data visualization are ultimately ones of design. I do believe that there are some hard and fast rules that should never be broken (e.g. always start the y-axis at 0 for bar charts of strictly positive values, don’t represent data with the same units on a secondary y-axis, never use a line chart for categorical data), however I also believe there are some that are more flexible, depending on what you want to accomplish, and your audience. Should you never, ever, use a pie chart? No. Some people are more comfortable with pie charts just from their familiarity with them. Is a bar chart a better choice in terms of representing the data? Yes. But that doesn’t mean there aren’t exceptions (just don’t make a 3D one).

The same individual that made the observation about dataviz not being taught also pointed out to me another factors that can influence design choices: what she called chart fatigue. Is the bar chart the best way to plot a single metric across a categorical variable? Almost always, yes. But show a room full of businesspeople bar chart after bar chart after bar chart and anyone can tell you that they’re all going to start to look the same, and interpretation of them is going to suffer as a result. Plus you’re probably going to lose the interest of your audience.

Practice makes perfect

In conclusion, I think that awareness of data visualization is only going to get better as companies (and the average consumer) become more “data savvy”. It is my sincere hope that people will give more and more emphasis, not only to the importance of visualization as a tool, but also to the design choices around it, and what constitutes good and bad depictions of data.

For now, just remember that data visualization is ultimately all about communicating and having your reader understand, not necessarily wowing them (though both together are not impossible). And sometimes, that means boring is better.

Visual Analytics of Every PS4 and XBox One Game by Install Size

Introduction

I have a startling confession to make: sometimes I just want things to be simple. Sometimes I just want things to be easier. All this talk about “Big Data”, predictive modeling, machine learning, and all those associated bits and pieces that go along with data science can be a bit mentally exhausting. I want to take a step back and work with a smaller dataset, something simple, something easy, something everyone can relate to – after all, that’s what this blog started out being about.

A while back, someone posted on Slashdot that the folks over at Finder.com had put together data sets of the install size of every PS4 and Xbox One game released to date. Being a a console owner myself – I’m a PS4 guy, but no fanboy or hardcore gamer by any means – I thought this would be a fun and rather unique data set to play around with, one that would fit well within the category of ‘everyday analytics’. So let’s take a look shall we?

Background

Very little background required here – the dataset comprises the title, release date, release type (major or indie), console (PS4 or Xbox One), and size in GiB of all games released as of September 10th, 2015. For this post we will ignore the time-related dimension and look only at the quantity of interest: install size.

Analysis

Okay, if I gave this data to your average Excel jockey what’s the first thing they’d do? A high level summary of the data broken apart by categorical variables and summarized by quantitative? You got it!
We can see that far more PS4 games have been released than Xbox (462 vs. 336) and the relative proportions are reversed for the former platform versus the latter as release type goes.

A small aside here on data visualization: it’s worth noting that the above is a good way to go for making a bar chart from a functional perspective. Since there are data labels and the y-axis metric is in the title, we can ditch the axis and maximize the data-ink ratio (well, data-pixel anyhow). I’ve also avoided using a stacked bar chart as interpretation of absolute values tends to suffer when not read from the same baseline. I’m okay with doing it for relative proportions though – as in the below, which further illustrates the difference in release type proportion between the two consoles:

Finally, how does the install size differ between the consoles and game types? If I’m an average analyst and just trying to get a grip on the data, I’d take an average to start:
We can see (unsurprisingly, if you know anything about console games) that major releases tend to be much larger in size than indie. Also in both cases, Xbox install sizes are larger on average: about 1.7x for indie titles and 1.25x for major.

Okay, that’s interesting. But if you’re like me, you’ll be thinking about how 99% of the phenomena in the universe are distributed by a power law or have some kind of non-Gaussian based distribution, and so averages are actually not always such a great way to summarize data. Is this the case for our install size data set?

Yes, it is. We can see here in this combination histogram / cumulative PDF (in the form of a Pareto chart) that the games follow a power law, with approximately 55 and 65 percent being < 5 GiB, for PS4 and Xbox games respectively

But is this entirely due to the indie games having small sizes? Might the major releases be centered around some average or median size?

No, we can see that even when broken apart by type of release the power-law like distribution for install sizes persists. I compared the averages to medians found them to be still be decent representations of central tendency and not too affected by outliers.

Finally we can look at the distribution of the install sizes by using another type of visualization suited for this task, the boxplot. While it is at least possible to jury-rig up a boxplot in Excel (see this excellent how-to over at Peltier Tech) Google Sheets doesn’t give us as much to work with, but I did my best (the data label is at the maximum point, and the value is the difference between the max and Q3):

The plots show that install sizes are generally greater for Xbox One vs. PS4, and that the difference (and skew) appears to be a bit more pronounced for indie games versus major releases, as we saw in the previous figures.

Okay, that’s all very interesting, but what about games that are available for both consoles? Are the install sizes generally the same or do they differ?

Difference in Install Size by Console Type
Because we’ve seen that the Xbox install sizes are generally larger than Playstation, here I take the PS4 size to be the baseline for games which appear on both (that is, differences are of the form XBOX Size – PS4 Size).

Of the 618 unique titles in the whole data set (798 titles if you double count across platform), 179 (~29%) were available on both – so roughly only a third of games are released for both major consoles.

Let’s take a look at the difference in install sizes – do those games which appear for both reflect what we saw earlier?

Yes, for both categories of game the majority are larger on Xbox than PS4 (none were the same size). Overall about 85% of the games were larger on Microsoft’s console (152/179).

Okay, but how much larger? Are we talking twice as large? Five times larger? Because the size of the games varies widely (especially between the release types) I opted to go for percentages here:

Unsurprisingly, on average indie games tend to have larger differences proportionally, because they’re generally much smaller in size than major releases. We can see they are nearly twice as large on Xbox vs. PS4 while major releases about 1 and a quarter. When games are larger on PS4, there’s not as big a disparity, and the pattern across release types is the same (though keep in mind the number of games here is a lot smaller than for the former).

Finally, just to ground this a bit more I thought I’d look at the top 10 games in each release type where the absolute differences are the largest. As I said before, here the difference is Xbox size minus PS4:

For major releases, the worst offender for being larger on PS4 is Batman: Arkham Night (~6.1 GiB difference) while on the opposite end, The Elder Scrolls Online has a ~19 GiB difference. Wow.

For indies, we can see the absolute difference is a lot smaller for those games bigger on PS4, with Octodad having the largest difference of ~1.4 GiB (56% of its PS4 size). Warframe is 19.6 GiB bigger on Xbox than PS4, or 503% larger (!!)

Finally, I’ve visualized all the data together for you so you can explore it yourself. Below is a bubble chart of the Xbox install size plotted against PS4, coloured by release type, where the size of each point represents the absolute value of the percentage difference between the platforms (with the PS4 size taken to be the baseline). So points above the diagonal are larger for Xbox than PS4, and points below the diagonal are larger for PS4 than Xbox. Also note that the scale is log-log. You can see that most of the major releases are pretty close to each other in size, as they nearly lie on the y=x line.

Conclusion

It’s been nice to get back into the swing of things and do a little simple data visualization, as well as play with a data set that falls into the ‘everyday analytics’ category.
And, as a result, we’ve learned:

  • XBox games generally tend to have larger install sizes than PS4 ones, even for the same title
  • Game install sizes follow a power law, just like most everything else in the universe (or maybe just 80% of it)
  • What the heck a GiB is
Until next time then, don’t fail to keep looking for the simple beauty in data all around you.

References & Resources

Complete List of Xbox One Install Sizes:
Complete List of PlayStation 4 Install Sizes:
Compiled data set (Google Sheets):
Excel Box and Whisker Diagrams (Box Plots) @ Peltier Tech:

Visualization and Analysis of Reddit’s “The Button” Data

Introduction

People are weird. And if there’s anything that’s greater collective proof of this fact than Reddit, you’d be hard pressed to find it.

I tend to put reddit in the same bucket as companies like Google, Amazon and Netflix, where they have enough money, or freedom, or both, to say something like “wouldn’t it be cool if….?” and then they do it simply because they can.

Enter “the button” (/r/thebutton), reddit’s great social experiment that appeared on April Fool’s Day of this year. An enticing blue rectangle with a timer that counts down from 60 to zero that’s reset when the button is pushed, with no explanation as to what happens when the time is allowed to run out. Sound familiar? The catch here being that it was an experience shared by anyone who visited the site, and each user also only got one press (though many made attempts to game the system, at least initially).

Finally, the timer reached zero, the last button press being at 2015-06-05 21:49:53.069000UTC, and the game (rather anti-climactically I might offer) ended.

What does this have to do with people being weird? Well, an entire mythology was built up around the button, amongst other things. Okay, maybe interesting is a better word. And maybe we’re just talking about your average redditor.

Either way, what interests me is that when the experiment ended, all the data were made available. So let’s have a look shall we?

Background

The dataset consists of simply four fields: 
press time, the date and time the button was pressed
flair, the flair the user was assigned given at what the timer was at when they pushed the button
css, the flair class given to the user
and lastly outage press, a Boolean indicator as to if the press occurred during a site outage.
The data span a time period from 2015-04-01 16:10:04.468000 to 2015-06-05 21:49:53.069000, with a total of 1,008,316 rows (unique presses).
I found there was css missing for some rows, and a lot of of “non presser” flair (users who were not eligible to press the button as their account was created after the event started). For these I used a “missing” value of -1 for the number of seconds remaining when the button was pushed; otherwise it could be stripped from the css field.

Analysis

With this data set, we’re looking at a pretty straightforward categorical time series.
Overall Activity in Time
First we can just look at the total number of button presses, regardless of what the clock said (when they occurred in the countdown) by plotting the raw number of presses per day:

Hmmm… you can see there is a massive spike at the beginning of the graph and there’s much, much fewer for the rest of the duration of the experiment. In fact, nearly 32% of all clicks occurred in the first day, and over half (51.3%) in the first two days. 
I think has something to do with both the initial interest in the experiment when it first was announced, and also with the fact that the higher the counter is kept at, the more people can press the button in the same time period (more on this later).
Perhaps a logarithmic graph for the y-axis would be more suitable?
That’s better. We can see the big drop-off in the first two days or so, and also that little blip around the 18th of May is more apparent. This is likely tied to one of several technical glitches which are noted in the button wiki,

For a more granular look, let’s do the hourly presses as well (with a log scale):

Cool. The spike on the 18th seems to be mainly around one hour with about a thousand presses, and we can see too that perhaps there’s some kind of periodic behavior in the data on an hourly basis? If we exclude some of the earlier data we can also go back to not using a log scale for the y-axis:

Let’s look more into the hours of the day when the button presses occur. We can create a simple bar plot of the count of button presses by hour overall:

You can see that the vast majority occurred around 5 PM and then there is a drop-off after that, with the lows being in the morning hours between about 7 and noon. Note that all the timestamps for the button pushes are in Universal Time. Unfortunately there is no geo data, but assuming most redditors who pushed the button are within the continental United States (a rather fair assumption) the high between 5-7 PM would be 11 AM to 1 PM (so, around your lunch hour at work).

But wait, that was just the overall sum of hours over the whole time period. Is there a daily pattern? What about by hour and day of week? Are most redditors pushing the button on the weekend or are they doing it at work (or during school)? We should look into this in more detail.

Hmm, nope! The majority of the clicks occurred Wednesday-Thursday night. But as we know from the previous graphs, the vast majority also occurred within the first two days, which happened to be a Wednesday and Thursday. So the figures above aren’t really that insightful, and perhaps it would make more sense to look at the trending in time across both day and hour? That would give us the figure as below:

As we saw before, there is a huge amount of clicks in the first few days (the first few hours even) so even with log scaling it’s hard to pick out a clear pattern. But most of the presses appear to be present in the bands after 15:00 and before 07:00. You can see the clicks around the outage on the 18th of May were in the same high period, around 18:00 and into the next day.

Maybe alternate colouring would help?

That’s better. Also if we exclude the flurry of activity in the first few days or so, we can drop the logarithmic scaling and see the other data in more detail:

To get a more normalized view, we can also look at the percentage of daily clicks per hour for each day, which yields a much more interesting view, and really shows the gap in the middle and the outage on the 18th:

Activity by Seconds Remaining
So far we’ve only looked at the button press activity by the counts in time. What about the time remaining for the presses? That’s what determined each individual reddit user’s flair, and was the basis for all the discussion around the button.

The reddit code granted flairs which were specific to the time remaining when the button was pushed.  For example, if there were 34 seconds remaining, then the css would be “34s”, so it was easy to strip these and convert into numeric data. There were also those that did not press the button who were given the “non presser” flair (6957 rows, ~0.69%), as well as a small number of entries missing flair (67, <0.01%), which I gave the placeholder value of -1.

The remaining flair classes served as a bucketing which functioned very much like a histogram:

Color Have they pressed? Can they press? Timer number when pressed
Grey/Gray N Y NA
Purple Y N 60.00 ~ 51.01
Blue Y N 51.00 ~ 41.01
Green Y N 41.00 ~ 31.01
Yellow Y N 31.00 ~ 21.01
Orange Y N 21.00 ~ 11.01
Red Y N 11.00 ~ 00.00
Silver/White N N NA

We can see this if we plot a histogram of the button presses by using the CSS class which gives the more granular seconds remaining, and use breaks the same as above:

We can see there is much greater proportion of those who pressed within 51-60s left, and there is falloff from there (power law). This is in line with what we saw in the time series graphs: the more the button was pressed, the more presses could occur in a given interval of time, and so we expect that most of those presses occurred during the peak activity at the beginning of the experiment (which we’ll soon examine).

What’s different from the documentation above from the button wiki is the “cheater” class, which was given to those who tried to game the system by doing things like disconnecting their internet and pressing the button multiple times (as far as I can tell). You can see that plotting a bar graph is similar to the above histogram with the difference being contained in the “cheater” class:

Furthermore, looking over the time period, how are the presses distributed in each class? What about in the cheater class? We can plot a more granular histogram:

Here we can more clearly see the exponential nature of the distribution, as well as little ‘bumps’ around the 10, 20, 30 and 45 second marks. Unfortunately this doesn’t tell us anything about the cheater class as it still has valid second values. So let’s do a boxplot by css class as well, showing both the classes (buckets) as well as their distributions:

Obviously each class has to fit into a certain range given their definition, but we can see some are more skewed than others (e.g. class for 51-60s is highly negatively skewed, whereas the class for 41-50 has median around 45). Also we can see that the majority of the cheater class is right near the 60 mark.

If we want to be fancier we can also plot the boxplot using just the points themselves and adding jitter:

This shows the skew of the distributions per class/bucket (focus around “round” times like 10, 30, 45s, etc.) as before, as well as how the vast majority of the cheater class appears to be at 59s mark.

Presses by seconds remaining and in time
Lastly we can combine the analyses above and look at how the quantity and proportion of button presses varies in time by the class and number of seconds remaining.

First we can look at the raw count of presses per css type per day as a line graph. Note again the scale on the y-axis is logarithmic:

This is a bit noisy, but we can see that the press-6 class (presses with 51-60s remaining) dominate at the beginning, then taper off toward the end. Presses in the 0-10 class did not appear until after April 15, then eventually overtook the quicker presses, as would have to be the case in order for the timer to run out. The cheater class starts very high with the press-6 class, then drops off significantly and continues to decrease. I would have like to break this out into small multiples for more clarity, but it’s not the easiest to do using ggplot.

Another way to look at it would be to look at the percent of presses by class per day. I’ve written previously about how stacked area graphs are not your friend, but in this case it’s actually not too bad (plus I wanted to learn how to do it in ggplot). If anything it shows the increase presses in the 51-60 range right after the outage on May 18, and the increase in the 0-10 range toward the end (green):

This is all very well and good, but let’s get more granular. We can easily visualize the data more granularly using heatmaps with the second values taken from the user flair to get a much more detailed picture. First we’ll look at a heatmap of this by hour over the time period:

Again, the scaling is logarithmic for the counts (here the fill colour). We can see some interesting patterns emerging, but it’s a little too sparse as there are a lot of hours without presses for a particular second value. Let’s really get granular and use all the data on the per second level!

On the left is the data for the whole period with a logartihmic scale, whereas the figure on the right excludes some of the earlier data and uses a linear scale. We can see the beginning peak activity in the upper lefthand corner, and then these interesting bands around the 5, 10, 20, 30, and 45 marks forming and gaining strength over time (particular toward the end). Interestingly in addition the resurgence in near-instantaneous presses after the outage around May 18, there was also a hotspot of presses around the 45s mark close to the end of April. Alternate colouring below:

Finally, we can divide by the number of presses per day and calculate the percent each number of seconds remaining made up over the time period. That gives the figures below:

Here the flurry of activity at the beginning continues to be prominent, but the bands also stand out a little more on a daily basis. We can also see how the proportion of clicks for the smaller number of seconds remaining continues to increase until finally the timer is allowed to run out.

Conclusion

The button experiment is over. In the end there was no momentous meaning to it all, no grand scheme or plan, no hatch exploding into the jungle, just an announcement that the thread would be archived. Again, somewhat anti-climactic.
But, it was an interesting experiment. This was an interesting data set, given the relationship between the amount of data that could exist in the same interval of time because of the nature of it. 
And I think it really says something about what the internet allows us to do (both in terms of creating something simply for the sake of it, and collecting and analyzing data), and also about people’s desire to find patterns and create meaning in things, no matter what they are. If you’d asked me, I never would have guessed religions would have sprung up around something as simple as pushing a button. But then again, religions have sprung up around stranger things.
You can read and discuss in the button aftermath thread, and if you want to have a go at it yourself, the code and data are below. Until next time I’ll just keep pressing on.

References & Resources

the button press data (from reddit’s github)

R code for plots

/r/thebutton

Data Visualization Fundamentals with Skittles

So I have a shocking confession to make: I love Skittles.

This post is not sponsored, endorsed, compensated, or paid for in any way, shape or form, by Skittles Candy. I’m not particular – I like other types of candy that are similar – you know, those ones that are chocolate covered in a hard shell, whether they be the kind where you eat the red ones last or not.

Anyhow, I got to thinking about how, abstractly, each individual candy can be viewed like a pixel of a different color. So you can make art using candy, just like artists make a mosaic. There’s lots of this on the internet you can already see: in fact, Skittles has done print advertising this way.

But…. each individual candy can also represent something else: a unit of measurement. I thought it would be cool to go through some data visualization fundamentals using the candy in this way. So let’s dive in.

Data Visualization using only 1 bag of Skittles

So, what would your average first grader do with a bag of Skittles if you asked them to sort it? Probably something like below, the physical equivalent of a bubble chart depicting the quantities of each colour by area, assuming each Skittle is approximately the same size.

A perhaps more useful way to do the same would be to organize each colour in rows, with each row a set number (like tally marks). Here it’s not only easy to see the relative proportions of the different colours in the bag, but also count them as each row and group is a set number (5 & 10, respectively). This is equivalent to a pictogram, with each Skittle representing, well, 1 Skittle:

It’s not a big stretch of the imagination to collapse those groups together into groups of a set height. So here we have a proportional bar chart, where the length of each bar represents the percentage of the bag that is each colour. Note that because I didn’t slice Skittles in half, the physical analogue is not exactly the same as what you’d put down on paper or in Excel (there is one additional unit for yellow and orange):

And, as I both often have to remind people of this rule, and also observe many people not following it, it is best practice to sort the bars in descending order for maximum clarity / comparative value (assuming there is not another more important ordering):

And, if we want to transform our proportional bar chart into one comparing absolute quantities, it is not a giant stretch of the imagination to break apart the different bars so they are only one ‘pixel’ high:

Here it’s much easier to get an idea of the absolute number of each colour in the bag, but harder to tally that numbers exactly – for that we’d need to add an axis or data labels.

More bags please

Okay, I have another shocking confession to make: I lied. I really like Skittles. So I actually bought a whole bunch of bags.

So let’s look at some more visualization fundamentals, where we required comparing not only across a categorical variable (colour) but also between groups.

Here is the equivalent to our first graph from before, only showing the different numbers of Skittles in each bag. You can see there’s actually a fair amount of variance; the smallest bag had 89 pieces of candy, whereas the largest had 110.

Now let’s make a bubble graph which not only compares the sizes between the different bags, but also their makeups by colour. The end result is actually closer to a collection of pie charts:

We can also group by colour only to see the overall makeup for the whole group of bags. Whereas orange dominated in the first bag we looked at, you can see here that orange and yellow are approximately at parity overall.

Now let’s look at the tally mark / pictograph method. Here each row represents a bag:

You can see there’s a fair bit of variance in the different colours. I also tried rearranging things so they result was less like a pictograph and more like a treemap:

Really the best way to compare would be a bar graph. Here’s a stacked area graph. I didn’t bother sorting by length, because at this point I was pretty tired of shuffling Skittles around:

To get a better idea of the different makeups of each bag by colour, we can break this out into a grouped bar graph, first by bag, then by colour:
And, of course, we can reverse the order if we want to more directly compare the colour makeups. The columns are in numerical order by bag. And just for fun, we’ll make this one a column chart:
There. That’s better! Clearly Bag 1 was an outlier as far as the number of purple went, and Bag 3 had a lot of yellow. 

Concluding Remark

I thought it’d be cool to mix things up a bit, and trying doing some data visualization using a physical medium. The end result ended up being something more like an exercise for an elementary school mathematics class (indeed, there are many examples of this online), but I think it still drives home some of the fundamental strengths and weakness of different visualization types, as well as showing how they can be depicted using different media. 
If you’re really interested, you can download the data yourself and slice and dice visualizations to your heart’s content. And I’m sure if you bought enough bags of Skittles you could learn something of a statistical nature about their manufacturing and packaging process – but perhaps that’s for a different day. Until then I’ll just enjoy good candy and data visualization.

Toronto Data Science Meetup – Machine Learning for Humans

A little while ago I spoke again at the Toronto Data Science Group, and gave a presentation I called “Machine Learning for Humans”:

I had originally intended to cover a wide variety of general “gotchas” around the practical applications of machine learning, however with half an hour there’s really only so much you can cover.

The talk ended up being more of an overview of binary classification, as well as some anecdotes around mistakes in using machine learning I’ve actually seen in the field, including:

  • Not doing any model evaluation at all
  • Doing model evaluation but without cross-validation
  • Not knowing what the cold start problem is and how to avoid it with a recommender system
All in all it was received very well despite being review for a lot of people in the room. As usual, I took away some learnings around presenting:
  • Always lowball for time (the presentation was rushed despite my blistering pace)
  • Never try to use fancy fonts in Powerpoint and expect them to carry over – it never works (copy paste as an image instead when you’ve got the final presentation)
Dan Thierl of Rubikloud gave a really informative and candid talk about what product management at a data science startup can look like. In particular, I was struck by his honesty around the challenges faced (both from technical standpoint and with clients), how quickly you have to move / pivot, and how some clients are just looking for simple solutions (Can you help us dashboard?) and are perhaps not at a level of maturity to want or fully utilize a data science solution.
All in all, another great meetup that prompted some really interesting discussion afterward. I look forward to the next one. I’ve added the presentation to the speaking section.

What a Gas! The Falling Price of Oil and Ontario Gasoline Prices

Introduction

In case you’ve been living under a rock, there’s been a lot of chatter in the financial world late about the price of oil going down. Way, way, down. So much so that the Bank of Canada cut interest rates. What crazy times are these we live in? I thought gas was only going to get more and more expensive until the end of time until everyone resorted to driving solar-powered cars.
I’m no economist (or commodities guy) and frankly, a lot of it seems like black magic and voodoo to me, but I thought it’d take a look at the data to see just how things have changed. Plus it’ll be an opportunity to muck about with some time series analysis in R.

Your average economists.

Background

Not much background is really needed other than what I mentioned in the introduction, but, well, we do need some data to work with.
Ontario Gas Prices for a number of cities (as well as the province averages) can be found at the Ontario Ministry of Energy. Unfortunately they’re year-by-year CSVs with a lot of different information row-wise (Weekly, Month, YTD).
No biggie, a simple bash script will take care of downloading all the data with wget and parsing and concatenating the data with other *nix tools:
# Download data
for i in $(seq 1990 2014)
do wget http://www.energy.gov.on.ca/fuelupload/ONTREG$i.csv
done

# Retain the header
head -n 2 ONTREG1990.csv | sed 1d > ONTREG_merged.csv

# Loop over the files and use sed to extract the relevant lines
for i in $(seq 1990 2014)
do
tail -n 15 ONTREG$i.csv | sed 13,15d | sed 's/./-01-'$i',/4' >> ONTREG_merged.csv
done
Great! Now we have all monthly Ontario Gas Price data from 1990-2014 inclusive in one big file.

The WTI data I got from The Federal Reserve Bank of St. Louis, and the forecasts from the US Energy Information Administration.

Analysis

First, a pretty plot of the data:

Wow, that really is a cliff, isn’t it? The average Ontario Gas Price hasn’t been as low as it was in Dec 2014 since the fall of 2010 (Sep 2010, $1.01) and the WTI not since about May of that year.

Now to the fun stuff. Let’s read the data into R and do some hacky time series analysis.

library(ts)
library(forecast)

# Read in the Ontario Gas Price Data
data <- read.csv(file="ONTREG_merged.csv", header=T, sep=",")

# Read in the WTI oil price data
WTI_data <- read.csv(file='DCOILWTICO.csv',header=F, col.names=c("Date", "Value"))

# Create a time series object for the WTI and Ontario Avg
WTI <- ts(data=WTI_data$Value, frequency=12, start=c(1990,1), end=c(2014,12))
ON <- ts(data=data$ON.Avg, frequency=12, start=c(1990,1), end=c(2014,12))

# Plot and compare
combined <- cbind(WTI, ON)
plot(combined)

We get a comparison plot of the data:

And we can look at the Ontario Gas Price as a function of the WTI. Linear is on the left, log-log on the right.

Next we build lm model objects and look at the diagnostics. I’ll spare the details, but I feel better about the log-log, so we’ll go with that.

# Create linear models (normal and log-log)
l1 <- lm(ON ~ WTI, data=combined)
l2 <- lm(log(ON) ~ log(WTI), data=combined)

# Compare relative performance
summary(l1)
summary(l2)
plot(l1)
plot(l2)

# Plot
plot(ON ~ WTI, data=combined, pch=16, cex=0.3)
abline(l1)
plot(log(ON) ~ log(WTI), data=combined, pch=16, cex=0.3)
abline(l2)

Lastly, we read in the forecast WTI data and use it to forecast the Ontario Gas price using our second model:
# Read in WTI forecast data
WTI_forecast <- read.csv(file="WTI_forecast.csv", header=F, sep=",", col.names=c("Date", "Value"))

# Forecast Ontario Gas Price
fit <- forecast(l2, newdata=data.frame(WTI=WTI_forecast$Value))

# Unlog
fit$mean <- exp(fit$mean)
fit$lower <- exp(fit$lower)
fit$upper <- exp(fit$upper)
fit$x <- exp(fit$x)

# Plot
plot(fit, ylab='Ontario Average Gas Price (cents/L)')

And there you have it! Putting the forecast data (blue line) and the WTI forecast back into our original graph, we can compare the two:

It’s that easy, right? That’s all there is to it?

Conclusion

Sadly, no. I’ve done some very quick work here and demonstrated some of the types of tools that are available in R, but “real” time-series analysis is a practice which requires much more care and nuance. 
For example, linear modeling assumes that variables are stationary (i.e. have constant mean and variance) and not auto-correlated, properties which are almost never true in the real world. Using methods such as above for non-stationary times series can result in what is known as “spurious regression” – finding relationships between variables which don’t really exist, even though the results have high R-squared and p-values.

In these cases, testing the stationary assumption and then massaging of the data (differencing & deseasonalization) is usually required beforehand to handle the problem, or other models which do not have as strict assumptions as linear regression are more appropriate. For more on this, see the excellent book “Forecasting: Principles and Practice” by Rob J. Hyndman and George Athana­sopou­los.

The take-away, as always, is that real life is complicated, and so analysis requires care. Predicting the future has never been, and never will be easy; either way, I just hope gas stays cheap.

References & Resources

code & data on github:

Fuel Price Data: Regular Unleaded Gasoline 2014 (Ontario Ministry of Energy):

http://www.energy.gov.on.ca/en/fuel-prices/fuel-price-data/?fuel=REG&yr=2015

Crude Oil Prices: West Texas Intermediate (WTI) (Federal Reserve Bank of St. Louis)

Short-Term Energy Outlook (U.S. Energy Information Administration)

Forecasting: Principles and Practice (OTexts)

The Mandelbrot Set in R

Introduction

I was having a conversation with someone I know about weather forecasts the other day and it went something like this:
“Yes, their forecasts are really not very good. They really need to work on improving them.”
“Well, yes, I think they’re okay, considering what they’re trying to do is impossible.”
The thought that a relatively simple set of equations can give rise to infinitely complex behavior is something which has intrigue and appeal beyond those mathematically minded and academic. Fractals became a part of pop culture, and the idea of chaos is not unknown to those outside mathematical research, as it was popularized by the term “The Butterfly Effect” (though unfortunately I can’t say the movie starring Ashton Kutcher was anything to write home about).
When I was finishing high school and starting my undergraduate degree I got very interested in fractals and chaos. What really spurred this interest was James Gleick’s Chaos, which outlined the history of the “mathematics of chaos” in a highly readable, narrative way. So much so that I later went on to take a pure mathematics course in chaos during my degree, and wrote a program in MATLAB for exploring the Mandelbrot set.
So I thought I’d give it a shot in R.

Background

Unquestionably the most well known fractal of all is the Mandelbrot fractal. Without going into laborious mathematical detail, the Mandelbrot set stems from using complex numbers, which are numbers of the form z = x + yi. The number i, the “imaginary unit”, is the special quantity such that i2 = –1. As such complex numbers are unique in that multiplying two of them together can result in much different behavior than numbers on the real number line: the magnitude of their product is not guaranteed to be greater than the terms multiplied, even for two quantities with real and imaginary parts with magnitude greater than one.
The Mandelbrot set the set of numbers produced by the following:
1. Initialize a complex set of numbers in the complex plane that are all z = 0.
2. Iterate the formula

such where c is a set of complex numbers filling the complex plane.
3. The Mandelbrot set is the set where z remains bounded for all n.
And mathematically it can be shown that if under iteration z is greater than 2 than c is not in the set.
In practice one keeps track of the number of iterations it takes z to diverge for a given c, and then colours the points accordingly to their “speed of divergence”.

Execution

In R this is straightforward code. Pick some parameters about the space (range and resolution for x and y, number of iterations, etc.) and then iterate. Naively, the code could be executed as below:
However, if you know anything about writing good code in R, you know that for loops are bad, and it’s better to rely upon the highly optimized underpinnings of the language by using array-based functions like which and lapply.  So a more efficient version of the same code is below:

We can then plunk these into functions where the different criteria for the rendered fractal are parameters:
Let’s compare the runtime between the two shall we? For the base settings, how does the naive version compare to the vectorized one? 
> compare_runtimes()
   user  system elapsed 
 37.713   0.064  37.773 
   user  system elapsed 
  0.482   0.094   0.576 
The results speak for themselves: a ~65x decrease in runtime by using R indexing functions instead of for loops! This definitely speaks to the importance of writing optimized code taking advantage of R’s vectorized functions.
You can tweak the different parameters for the function to return different parts of the complex space. Below are some pretty example plots of what the script can do with different parameters:

Conclusion

What does all this have to do with my conversation about the weather? Why, everything, of course! It’s where chaos theory came from. Enjoy the pretty pictures. It was fun getting R to churn out some nice fractals and good to take a trip down memory lane.

References and Resources

The Mandelbot Set: 
code on github:

Toronto Cats and Dogs II – Top 25 Names of 2014

I was quite surprised by the relative popularity of my previous analysis of the data for Licensed Cats & Dogs in Toronto for 2011, given how simple it was to put together.

I was browsing the Open Data Portal recently and noticed that there was a new data set for pets: the top 25 names for both dogs and cats. I thought this could lend itself to some quick, easy visualization and be a neat little addition to the previous post.

First we simply visualize the raw counts of the Top 25 names against each other. Interestingly, the top 2 names for both dogs and cats are apparently the same: Charlie and Max.

Next let’s take a look at the distribution of these top 25 names for each type of pet by how long they are, which just involves calculating the name length and then pooling the counts:

You can see that, proportionally the top dog names are a bit shorter (distribution is positively / right-skewed) compared to the cat names (slightly negatively / left skewed). Also note both are centered around names of length 5, and the one cat name of length 8 (Princess).

Looking at the dog names, do you notice something interesting about them? A particular feature present in nearly all? I did. Nearly every one of the top 25 dog names ends in a vowel. We can see this by visualizing the proportion of the counts for each type of pet by whether the name ends in a vowel or consonant:

Which to me, seems to indicate that more dogs tend to have “cutesy” names, usually ending in ‘y’, than cats.

Fun stuff, but one thing really bothers me… no “Fido” or “Boots”? I guess some once popular names have gone to the dogs.

References & Resources

Licensed Dog and Cat Names (Toronto Open Data)

Twitter Pop-up Analytics

Introduction

So I’ve been thinking a lot lately. Well, that’s always true. I should say, I’ve been thinking a lot lately about the blog. When I started this blog I was very much into the whole quantified self thing, because it was new to me, I liked the data collection and analysis aspect, and I had a lot of time to play around with these little side projects.
When I started the blog I called it “everyday analytics” because that’s what I saw it always being; analysis of data on topics that were part of everyday life, the ordinary viewed under the analytical lens, things that everyone can relate to. You can see this in my original about page for the blog which has remained the same since inception.
I was thinking a lot lately about how as my interest in data analysis, visualization and analytics has matured, and so that’s not really the case so much anymore. The content of everyday analytics has become a lot less everyday. Analyzing the relative nutritional value of different items on the McDonald’s menu (yeesh, looking back now those graphs are pretty bad) is very much something to which most everyone could relate. 2-D Histograms in R? PCA and K-means clustering? Not so much.
So along this line of thinking, for this reason, I thought it’s high time to get back into the original spirit of the site when it was started. So I thought I’d do some quick quantified-self type analysis, about something everyone can relate to, nothing fancy. 
Let’s look at my Twitter feed.

Background

It wasn’t always easy to get data out of Twitter. If you look back at how Twitter’s API has changed over the years, there has been considerable uproar about the restrictions they’ve made in updates, however they’re entitled to do so as they do hold the keys to the kingdom after all (it is their product). In fact, I thought it’d be a easiest to do this analysis just using the twitteR package, but it appears to be broken since Twitter has made said updates to their API.

Luckily I am not a developer. My data needs are simple for some ad hoc analysis. All I need is the data pulled and I am ready to go. Twitter now makes this easy now for anyone to do, just go to your settings page:

And then select the ‘Download archive’ button under ‘Your Twitter Archive’ (here it is a prompt to resend mine, as I took the screenshot after):

And boom! A CSV of all your tweets is in your inbox ready for analysis. After all this talk about working with “Big Data” and trawling through large datasets, it’s nice to take a breather a work with something small and simple.

Analysis

So, as I said, nothing fancy here, just wrote some intentionally hacky R code to do some “pop-up” analytics given Twitter’s output CSV. Why did I do it this way, which results in 1990ish looking graphs, instead of in Excel and making it all pretty? Why, for you, of course. Reproducibility. You can take my same R code and run it on your twitter archive (which is probably a lot larger and more interesting than mine) and get the same graphs.
The data set comprises 328 tweets sent by myself between 2012-06-03 and 2014-10-02. The fields I examined were the datetime field (time parting analysis), the tweet source and the text / content.
Time Parting
First let’s look at the time trending of my tweeting behaviour:
We can see there is some kind of periodicity, with peaks and valleys in how many tweets I send. The sharp decline near the end is because there are only 2 days of data for October. Also, compared to your average Twitter user, I’d say I don’t tweet alot, generally only once every two days or so on average:
> summary(as.vector(monthly))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1.00    8.00   12.00   11.31   15.00   21.00 

Let’s take a look and see if there is any rhyme or reason to these peaks and valleys:

Looking at the total counts per month, it looks like I’ve tweeted less often in March, July and December for whatever reason (for all of this, pardon my eyeballing..)

What about by day of week?

Look like I’ve tweeted quite a bit more on Tuesday, and markedly less on the weekend. Now, how does that look over the course of the day;

My peak tweeting time seems to be around 4 PM. Apparently I have sent tweets even in the wee hours of the morning – this was a surprise to me. I took a stab at making a heatmap, but it was quite sparse; however the 4-6 PM peak does persist across the days of the week.

Tweets by Source
Okay, that was interesting. Where am I tweeting from?

Look like the majority of my tweets are actually sent from the desktop site, followed by my phone, and then sharing on sites. I attribute this to the fact that I mainly use twitter to share articles, which isn’t easy to do on my smartphone.

Content Analysis
Ah, now on to the interesting stuff! What’s actually in those tweets?
First let’s look at the length of my tweets in a simple histogram:
Looks like generally my tweets are above 70 characters or so, with a large peak close to the absolute limit of 160 characters. 
Okay, but what I am actually tweeting about? Using the very awesome tm package it’s easy to do some simple text mining and pull out both top frequent terms, as well as hashtags.
So apparently I tweet a lot about data, analysis, Toronto and visualization. To anyone who’s read the blog this shouldn’t be overly surprisingly. Also you can see I pass along articles and interact with others as “via” and “thanks” are in there too. Too bad about that garbage ampersand.
Overwhelmingly the top hashtag I use is #dataviz, followed of course by #rstats. Again, for anyone who knows me (or has seen one of my talks) this should not come as a surprise. You can also see my use of Toronto Open Data in the #opendata and #dataeh hashtags.

Conclusion

That’s all for now. As I said, this was just a fun exercise to write some quick, easy R code to do some simple personal analytics on a small dataset. On the plus side the code is generalized, so I invite you to take it and look at your own twitter archive.
Or, you could pull all of someone else’s tweets, but that would, of course, require a little more work.

References

code at github
Twitter Help Center: Downloading Your Archive
The R Text Mining ™ package at CRAN
twitteR package at CRN

5 Ways to Do 2D Histograms in R

Introduction

Lately I was trying to put together some 2D histograms in R and found that there are many ways to do it, with directions on how to do so scattered across the internet in blogs, forums and of course, Stackoverflow.

As such I thought I’d give each a go and also put all of them together here for easy reference while also highlighting their difference.

For those not “in the know” a 2D histogram is an extensions of the regular old histogram, showing the distribution of values in a data set across the range of two quantitative variables. It can be considered a special case of the heat map, where the intensity values are just the count of observations in the data set within a particular area of the 2D space (bucket or bin).

So, quickly, here are 5 ways to make 2D histograms in R, plus one additional figure which is pretty neat.

First and foremost I get the palette looking all pretty using RColorBrewer, and then chuck some normally distributed data into a data frame (because I’m lazy). Also one scatterplot to justify the use of histograms.

# Color housekeeping
library(RColorBrewer)
rf <- colorRampPalette(rev(brewer.pal(11,'Spectral')))
r <- rf(32)

# Create normally distributed data for plotting
x <- rnorm(mean=1.5, 5000)
y <- rnorm(mean=1.6, 5000)
df <- data.frame(x,y)

# Plot
plot(df, pch=16, col='black', cex=0.5)

Option 1: hexbin

The hexbin package slices the space into 2D hexagons and then counts the number of points in each hexagon. The nice thing about hexbin is that it provides a legend for you, which adding manually in R is always a pain. The default invocation provides a pretty sparse looking monochrome figure. Adding the colramp parameter with a suitable vector produced from colorRampPalette makes things nicer. The legend placement is a bit strange – I adjusted it after the fact though you just as well do so in the R code.
##### OPTION 1: hexbin from package 'hexbin' #######
library(hexbin)
# Create hexbin object and plot
h <- hexbin(df)
plot(h)
plot(h, colramp=rf)

Using the hexbinplot function provides greater flexibility, allowing specification of endpoints for the bin counting, and also allowing the provision of a transformation function. Here I did log scaling. Also it appears to handle the legend placement better; no adjustment was required for these figures.

# hexbinplot function allows greater flexibility
hexbinplot(y~x, data=df, colramp=rf)
# Setting max and mins
hexbinplot(y~x, data=df, colramp=rf, mincnt=2, maxcnt=60)

# Scaling of legend - must provide both trans and inv functions
hexbinplot(y~x, data=df, colramp=rf, trans=log, inv=exp)

Option 2: hist2d

Another simple way to get a quick 2D histogram is to use the hist2d function from the gplots package. Again, the default invocation leaves a lot to be desired:
##### OPTION 2: hist2d from package 'gplots' #######
library(gplots)

# Default call
h2 <- hist2d(df)
Setting the colors and adjusting the bin sizing coarser yields a more desirable result. We can also scale so that the intensity is logarithmic as before.
# Coarser binsizing and add colouring
h2 <- hist2d(df, nbins=25, col=r)

# Scaling with log as before
h2 <- hist2d(df, nbins=25, col=r, FUN=function(x) log(length(x)))

Option 3: stat_2dbin from ggplot

And of course, where would a good R article be without reference to the ggplot way to do things? Here we can use the stat_bin2d function, either added to a ggplot object or as a type of geometry in the call to qplot.
##### OPTION 3: stat_bin2d from package 'ggplot' #######
library(ggplot2)

# Default call (as object)
p <- ggplot(df, aes(x,y))
h3 <- p + stat_bin2d()
h3

# Default call (using qplot)
qplot(x,y,data=df, geom='bin2d')
Again, we probably want to adjust the bin sizes to a desired number, and also ensure that ggplot uses our colours that we created before. The latter is done by adding the scale_fill_gradientn function with our colour vector as the colours argument. Log scaling is also easy to add using the trans parameter.
# Add colouring and change bins
h3 <- p + stat_bin2d(bins=25) + scale_fill_gradientn(colours=r)
h3

# Log scaling
h3 <- p + stat_bin2d(bins=25) + scale_fill_gradientn(colours=r, trans="log")
h3

Option 4: kde2d

Option #4 is to do kernel density estimation using kde2d from the MASS library. Here we are actually starting to stray from discrete bucketing of histograms to true density estimation, as this function does interpolation.
The default invocation uses n = 25 which is actually what we’ve been going with in this case. You can then plot the output using image().

Setting n higher does interpolation and we are into the realm of kernel density estimation, as you can set your “bin size” lower than how your data actually appear. Hadley Wickham notes that in R there are over 20 packages [PDF] with which to do density estimation so we’ll keep that to a separate discussion.

##### OPTION 4: kde2d from package 'MASS' #######
# Not a true heatmap as interpolated (kernel density estimation)
library(MASS)

# Default call
k <- kde2d(df$x, df$y)
image(k, col=r)

# Adjust binning (interpolate - can be computationally intensive for large datasets)
k <- kde2d(df$x, df$y, n=200)
image(k, col=r)

Option 5: The Hard Way

Lastly, an intrepid R user was nice enough to show on Stackoverflow how do it “the hard way” using base packages.
##### OPTION 5: The Hard Way (DIY) #######
# http://stackoverflow.com/questions/18089752/r-generate-2d-histogram-from-raw-data
nbins <- 25
x.bin <- seq(floor(min(df[,1])), ceiling(max(df[,1])), length=nbins)
y.bin <- seq(floor(min(df[,2])), ceiling(max(df[,2])), length=nbins)

freq <- as.data.frame(table(findInterval(df[,1], x.bin),findInterval(df[,2], y.bin)))
freq[,1] <- as.numeric(freq[,1])
freq[,2] <- as.numeric(freq[,2])

freq2D <- diag(nbins)*0
freq2D[cbind(freq[,1], freq[,2])] <- freq[,3]

# Normal
image(x.bin, y.bin, freq2D, col=r)

# Log
image(x.bin, y.bin, log(freq2D), col=r)
Not the way I would do it, given all the other options available, however if you want things “just so” maybe it’s for you.

Bonus Figure

Lastly I thought I would include this one very cool figure from Computational Actuarial Science with R which is not often seen, which includes both a 2D histogram with regular 1D histograms bordering it showing the density across each dimension.
##### Addendum: 2D Histogram + 1D on sides (from Computational ActSci w R) #######
#http://books.google.ca/books?id=YWcLBAAAQBAJ&pg=PA60&lpg=PA60&dq=kde2d+log&source=bl&ots=7AB-RAoMqY&sig=gFaHSoQCoGMXrR9BTaLOdCs198U&hl=en&sa=X&ei=8mQDVPqtMsi4ggSRnILQDw&redir_esc=y#v=onepage&q=kde2d%20log&f=false

h1 <- hist(df$x, breaks=25, plot=F)
h2 <- hist(df$y, breaks=25, plot=F)
top <- max(h1$counts, h2$counts)
k <- kde2d(df$x, df$y, n=25)

# margins
oldpar <- par()
par(mar=c(3,3,1,1))
layout(matrix(c(2,0,1,3),2,2,byrow=T),c(3,1), c(1,3))
image(k, col=r) #plot the image
par(mar=c(0,2,1,0))
barplot(h1$counts, axes=F, ylim=c(0, top), space=0, col='red')
par(mar=c(2,0,0.5,1))
barplot(h2$counts, axes=F, xlim=c(0, top), space=0, col='red', horiz=T)

Conclusion

So there you have it! 5 ways to create 2D histograms in R, plus some additional code to create a really snappy looking figure which incorporates the regular variety. I leave it to you to write (or find) some good code for creating legends for those functions which do not include them. Hopefully other R users will find this a helpful reference.

References

code on github
R generate 2D histogram from raw data (Stackoverflow)
Computational Actuarial Science with R (Google Books)
Wickham, Hadley. Density Estimation in R [PDF]