<h2>Abstract</h2>
<p>We perform a social experiment to investigate, if zombie related twitter posts can used as a reliable indicator for an early warning system. We show how such a system can be set up almost out-of-the-box using R - a free software environment for statistical computing and graphics. <strong>Warning</strong>: This blog entry contains toxic doses of Danish irony and sarcasm as well as disturbing graphs.</p>
<p><br><a rel="license" href="http://creativecommons.org/licenses/by-sa/4.0/"></a> This work is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by-sa/4.0/">Creative Commons Attribution-ShareAlike 4.0 International License</a>. The markdown+Rknitr source code of this blog is available under a <a href="https://www.gnu.org/licenses/gpl-3.0.html">GNU General Public License (GPL v3)</a> license from <a href="https://github.com/hoehleatsu/hoehleatsu.github.io"></a>.</p>
<h2>Introduction</h2>
<p>Proposing statistical methods is only mediocre fun if nobody applies them. As an act of desperation the prudent statistician has been forced to provide R packages supplemented with a CRAN, github, useR! or word-of-mouth advertising strategy. To underpin efforts, a reproducibility-crisis has been announced in order to scare decent comma-separated scientist <a href="https://www.washingtonpost.com/news/wonk/wp/2016/08/26/an-alarming-number-of-scientific-papers-contain-excel-errors/">from using Excel</a>. Social media marketing strategies of your R package include hashtag <code>#rstats</code> twitter announcements, possibly enhanced by a picture or animation showing your package at its best:</p>
<blockquote class="twitter-tweet">
<p dir="ltr" lang="en">
Introducing gganimate: <a href="https://twitter.com/hashtag/rstats?src=hash">#rstats</a> package for adding animation to any ggplot2 figure <a href="https://t.co/UBWKHmIc0e">https://t.co/UBWKHmIc0e</a> <a href="https://t.co/oQhQaYBqOj">pic.twitter.com/oQhQaYBqOj</a>
</p>
— David Robinson (@drob) <a href="https://twitter.com/drob/status/694274942813102080">February 1, 2016</a>
</blockquote>
<p>
</p><p>Unfortunately, little experience with the interactive aspect of this statistical software marketing strategy appears to be available. In order to fill this scientific advertising gap, this blog post constitutes an advertisement for the <strong>out-of-the-box-functionality</strong> of the <code>surveillance</code> package hidden as social experiment. It shows shows what you can do with R when combining a couple of packages, wrangle the data, cleverly visualize the results and then team up with the fantastic R community.</p>
<h2>The Setup: Detecting a Zombie Attack</h2>
<p>As previously explained in an <a href="http://user2015.math.aau.dk/lightning_talks">useR! 2015 lightning talk</a>, Max Brooks' <a href="https://en.wikipedia.org/wiki/The_Zombie_Survival_Guide">Zombie Survival Guide</a> is very concerned about the <strong>early warning</strong> of Zombie outbreaks.</p>
<a href="http://staff.math.su.se/hoehle/software/surveillance/hoehle-userR2015-web.pdf"></a>
<p><br></p>
<p>However, despite of extensive research and recommendations, no reliable service appears available for the early detection of such upcoming events. Twitter, on the other hand, has become the media darling to stay informed about news as they unfold. Hence, continuous monitoring of hashtags like <code>#zombie</code> or <code>#zombieattack</code> appears an essential component of your zombie survival strategy.</p>
<h1>Tight Clothes, Short Hair and R</h1>
<p>Extending the recommendations of the Zombie Survival guide we provide an out-of-the-box (OOTB) monitoring system by using the <code>rtweet</code> R package to obtain all individual tweets containing the hashtags <code>#zombie</code> or <code>#zombieattack</code>.</p>
<div><pre><code>the_query &lt;-<span> &quot;#zombieattack OR #zombie&quot;</span>
geocode &lt;-<span> &quot;&quot;</span>  <span>#To limit the seach to berlin &amp; surroundings: geocode &lt;- &quot;52.520583,13.402765,25km&quot;</span>
<span>#Converted query string which works for storing as file</span>
safe_query &lt;-<span> </span>stringr::<span>str_replace_all</span>(the_query, <span>&quot;[^[:alnum:]]&quot;</span>, <span>&quot;X&quot;</span>)</code></pre></div>
<p>In particular, the <a href="https://github.com/mkearney/rtweet">README</a> of the <code>rtweet</code> package provides helpful information on how to create a twitter app to automatically search tweets using the twitter API. One annoyance of the twitter REST API is that only the tweets of the past 7 days are kept in the index. Hence, your time series are going to be short unless you accumulate data over several queries spread over a time period. Instead of using a fancy database setup for this data collection, we provide a simple R solution based on <code>dplyr</code> and <code>saveRDS</code> - see the underlying R code of this post by clicking on the github logo in the license statement of this post. Basically,</p>
<ul><li>all tweets fulfilling the above hashtag search queries are extracted</li>
<li>each tweet is extended with a time stamp of the query-time</li>
<li>the entire result of each query us stored into a separate RDS-files using <code>saveRDS</code></li>
</ul><p>In a next step, all stored queries are loaded from the RDS files and put together. Subsequently, only the newest time stamped entry about each tweet is kept - this ensures that the re-tweeted counts are up-to-date and no post is counted twice. All these data wrangling operations are easily conducted using <code>dplyr</code>. Of course a full database solution would have been more elegant, but R does the job just as well as long it's not millions of queries. No matter the data backend, at the end of this pipeline we have a database of tweets.</p>
<div><pre><code><span>#Read the tweet database</span>
tw &lt;-<span> </span><span>readRDS</span>(<span>file=</span><span>paste0</span>(filePath,<span>&quot;Tweets-Database-&quot;</span>,safe_query,<span>&quot;-&quot;</span>,<span>&quot;2016-09-25&quot;</span>,<span>&quot;.RDS&quot;</span>))
<span>options</span>(<span>width=</span><span>300</span>,<span>tibble.width =</span> <span>Inf</span>)
tw %>%<span> </span><span>select</span>(created_at, retweet_count,screen_name,text,hashtags,query_at)</code></pre></div>
<pre><code>## # A tibble: 10,974 × 6
##             created_at retweet_count    screen_name                                                                                                                                          text  hashtags            query_at
##                 &lt;dttm>         &lt;int>          &lt;chr>                                                                                                                                         &lt;chr>    &lt;list>              &lt;dttm>
## 1  2016-09-25 10:26:28             0       Lovebian                                               The latest #Zombie Nation! https://t.co/8ZkOFSZH2v Thanks to @NJTVNews @MaxfireXSA @Xtopgun901X &lt;chr [1]> 2016-09-25 10:30:44
## 2  2016-09-25 10:25:49             2  MilesssAwaaay RT @Shaaooun: I'm gonna turn to a zombie soon! xdxdxdxd #AlmostSurvived #204Days #ITried #Zombie #StuckInMyRoom #Haha\n\n#MediaDoomsDay #Kame &lt;chr [7]> 2016-09-25 10:30:44
## 3  2016-09-25 10:21:10             6 catZzinthecity          RT @ZombieEventsUK: 7 reasons #TheGirlWithAllTheGifts is the best #zombie movie in years https://t.co/MB82ssxss2 via @MetroUK #Metro &lt;chr [3]> 2016-09-25 10:30:44
## 4  2016-09-25 10:19:41             0  CoolStuff2Get                             Think Geek Zombie Plush Slippers https://t.co/0em920WCMh #Zombie #Slippers #MyFeetAreCold https://t.co/iCEkPBykCa &lt;chr [3]> 2016-09-25 10:30:44
## 5  2016-09-25 10:19:41             4  TwitchersNews    RT @zOOkerx: Nur der frhe Vogel fngt den #zombie also schaut gemtlich rein bei @booty_pax! Now live #dayz on #twitch \n\nhttps://t.co/OIk6 &lt;chr [3]> 2016-09-25 10:30:44
## 6  2016-09-25 10:17:45             0 ZombieExaminer     Washington mall shooting suspect Arcan Cetin was '#Zombie-like' during arrest - USA TODAY https://t.co/itoDXG3L8T https://t.co/q2mURi24DB &lt;chr [1]> 2016-09-25 10:30:44
## 7  2016-09-25 10:17:44             4       SpawnRTs    RT @zOOkerx: Nur der frhe Vogel fngt den #zombie also schaut gemtlich rein bei @booty_pax! Now live #dayz on #twitch \n\nhttps://t.co/OIk6 &lt;chr [3]> 2016-09-25 10:30:44
## 8  2016-09-25 10:17:23             0   BennyPrabowo                   bad miku - bad oni-chan... no mercy\n.\n.\n.\n.\n#left4dead #games #hatsunemiku #fps #zombie #witch https://t.co/YP0nRDFFj7 &lt;chr [6]> 2016-09-25 10:30:44
## 9  2016-09-25 10:12:53            62   Nblackthorne  RT @PennilessScribe: He would end her pain, but he could no longer live in a world that demanded such sacrifice. #zombie #apocalypse\nhttps: &lt;chr [2]> 2016-09-25 10:30:44
## 10 2016-09-25 10:06:46             0   mthvillaalva                                                             Pak ganern!!! Kakatapos ko lang kumain ng dugo! \n#Zombie https://t.co/Zyd0btVJH4 &lt;chr [1]> 2016-09-25 10:30:44
## # ... with 10,964 more rows</code></pre>
<h3>OOTB Zombie Surveillance</h3>
<p>We are now ready to prospectively detect changes using the <code>surveillance</code> R package <span>(Salmon, Schumacher, and Höhle 2016)</span>.</p>
<div><pre><code><span>library</span>(<span>&quot;surveillance&quot;</span>)</code></pre></div>
<p>We shall initially focus on the <code>#zombie</code> series as it contains more counts. The first step is to convert the <code>data.frame</code> of individual tweets into a time series of daily counts.</p>
<div><pre><code><span>#' Function to convert data.frame to queries. For convenience we store the time series</span>
<span>#' and the data.frame jointly as a list. This allows for easy manipulations later on</span>
<span>#' as we see data.frame and time series to be a joint package.</span>
<span>#'</span>
<span>#' @param tw data.frame containing the linelist of tweets.</span>
<span>#' @param the_query_subset String containing a regexp to restrict the hashtags</span>
<span>#' @return List containing sts object as well as the original data frame.</span>
<span>#'</span>
df_2_timeseries &lt;-<span> </span>function(tw, the_query_subset) {
  tw_subset &lt;-<span> </span>tw %>%<span> </span><span>filter</span>(<span>grepl</span>(<span>gsub</span>(<span>&quot;#&quot;</span>,<span>&quot;&quot;</span>,the_query_subset),hashtags))

  <span>#Aggregate data per day and convert times series to sts object</span>
  ts &lt;-<span> </span>surveillance::<span>linelist2sts</span>(<span>as.data.frame</span>(tw_subset), <span>dateCol=</span><span>&quot;created_at_Date&quot;</span>, <span>aggregate.by=</span><span>&quot;1 day&quot;</span>)
  <span>#Drop first day with observations, due to the moving window of the twitter index, this count is incomplete</span>
  ts &lt;-<span> </span>ts[-<span>1</span>,]

  <span>return</span>(<span>list</span>(<span>tw=</span>tw_subset,<span>ts=</span>ts, <span>the_query_subset=</span>the_query_subset))
}

zombie &lt;-<span> </span><span>df_2_timeseries</span>(tw, <span>the_query_subset =</span> <span>&quot;#zombie&quot;</span>)</code></pre></div>
<p>It's easy to visualize the resulting time series using the plotting functionality of the surveillance package.</p>
<p>We see that the counts on the last day are incomplete. This is because the query was performed at 10:30 CEST and not at midnight. We therefore adjust counts on the last day based on simple inverse probability weighting. This just means that we scale up the counts by the inverse of the fraction the query-hour (10:30 CEST) makes up of 24h (see github code for details). This relies on the assumption that queries are evenly distributed over the day.</p>
<p>We are now ready to apply a surveillance algorithm to the pre-processed time series. We shall pick the so called C1 version of the EARS algorithm documented in <span>Hutwagner et al. (2003)</span> or <span>Fricker, Hegler, and Dunfee (2008)</span>. For a monitored time point <span>\(s\)</span> (here: a particular day, say 2016-09-23), this simple algorithm takes the previous seven observations before <span>\(s\)</span> in order to compute the mean and standard deviation, i.e. <span>\[
\begin{align*}
\bar{y}_s &amp;= \frac{1}{7} \sum_{t=s-8}^{s-1} y_t, \\
\operatorname{sd}_s &amp;= \frac{1}{7-1} \sum_{t=s-8}^{s-1} (y_t - \bar{y}_s)^2
\end{align*}
\]</span> The algorithm then computes the z-statistic <span>\(\operatorname{C1}_s = (y_s - \bar{y}_s)/\operatorname{sd}_s\)</span> for each time point to monitor. Once the value of this statistic is above 3 an alarm is flagged. This means that we assume that the previous 7 observations are what is to be expected when no unusual activity is going on. One can interpret the statistic as a transformation to (standard) normality: once the current observation is too extreme under this model an alarm is sounded. Such normal-approximations are justified given the large number of daily counts in the zombie series we consider, but does not take secular trends or day of the week effects into account. Note that the calculations can also be reversed in order to determine how large the number of observations need to be in order to generate an alarm.</p>
<p>We now apply the EARS C1 monitoring procedure to the zombie time series starting at the 8th day of the time series. It is important to realize that the result of monitoring a time point in the graphic is obtained by only <strong>looking into the past</strong>. Hence, the relevant time point to consider today is if an alarm would have occurred 2016-09-25. We also show the other time points to see, if we could have detected potential alarms earlier.</p>
<div><pre><code>zombie[[<span>&quot;sts&quot;</span>]] &lt;-<span> </span><span>earsC</span>(zombie$ts, <span>control=</span><span>list</span>(<span>range =</span> <span>8</span>:<span>nrow</span>(zombie$ts),
                         <span>method =</span> <span>&quot;C1&quot;</span>, <span>alpha =</span> <span>1</span>-<span>pnorm</span>(<span>3</span>)))</code></pre></div>
<p>What a relief! No suspicious zombie activity appears to be ongoing. Actually, it would have taken 511 tweets before we would have raised an alarm on 2016-09-25. This is quite a number.</p>
<p>As an additional sensitivity analysis we redo the analyses for the <code>#zombieattack</code> hashtag. Here the use of the normal approximation in the computation of the alerts is more questionable. Still, we can get a time series of counts together with the alarm limits.</p>
<p>Also no indication of zombie activity. The number of additional tweets needed before alarm in this case is: 21. Altogether, it looks safe out there...</p>
<h2>Summary</h2>
<p>R provides ideal functionality to quickly extract and monitor twitter time series. Combining with statistical process control methods allows you to prospectively monitor the use of hashtags. Twitter has released a dedicated package for this purpose, however, in case of low count time series it is better to use count-time series monitoring devices as implemented in the <code>surveillance</code> package. <span>Salmon, Schumacher, and Höhle (2016)</span> contains further details on how to proceed in this case.</p>
<p>The important question although remains: Does this really work in practice? Can you sleep tight, while your R zombie monitor scans twitter? Here is where the <strong>social experiment starts</strong>: Please help answer this question by retweeting the post below to create a drill alarm situation. More than 511 (!) and 21 tweets, respectively, are needed before an alarm will sound.</p>
(placeholder tweet, this will change in a couple of minutes!!)
<blockquote class="twitter-tweet">
<p dir="ltr" lang="en">
Video recording, slides &amp; R code of our <a href="https://twitter.com/ISDS"><span>(<span><strong>???</strong></span>)</span></a> MV Time Series webinar now available at <a href="https://t.co/XVtLrjbJKZ">https://t.co/XVtLrjbJKZ</a> <a href="https://twitter.com/hashtag/biosurveillance?src=hash">#biosurveillance</a> <a href="https://twitter.com/hashtag/rstats?src=hash">#rstats</a>
</p>
— Michael Höhle (@m_hoehle) <a href="https://twitter.com/m_hoehle/status/778712624426995712">21. September 2016</a>
</blockquote>
<p>I will continuously update the graphs in this post to see how our efforts are reflected in the time series of tweets containing the <code>#zombieattack</code> and <code>#zombie</code> hashtags. Thanks for your help!</p>
<p>
</p>
<p>
</p><h1>References</h1>
<div>
<div>
<p>Fricker, R. D., B. L. Hegler, and D. A. Dunfee. 2008. “Comparing syndromic surveillance detection methods: EARS’ versus a CUSUM-based methodology.” Stat Med 27 (17): 3407–29.</p>
</div>
<div>
<p>Hutwagner, L., W. Thompson, G. M. Seeman, and T. Treadwell. 2003. “The bioterrorism preparedness and response Early Aberration Reporting System (EARS).” J Urban Health 80 (2 Suppl 1): 89–96.</p>
</div>
<div>
<p>Salmon, M., D. Schumacher, and M. Höhle. 2016. “Monitoring Count Time Series in R: Aberration Detection in Public Health Surveillance.” Journal of Statistical Software 70 (10). doi:<a href="https://doi.org/10.18637/jss.v070.i10">10.18637/jss.v070.i10</a>.</p>
</div></div>