<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/combined-outlier-detection-with-dplyr-and-ruler/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://www.questionflow.org/2017/12/26/combined-outlier-detection-with-dplyr-and-ruler/"> QuestionFlow </a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p> Overview of simple outlier detection methods with their combination using dplyr and ruler packages.</p>
<div>
<h1>Prologue</h1>
<p>During the process of data analysis one of the most crucial steps is to identify and account for <a rel="nofollow" href="https://en.wikipedia.org/wiki/Outlier" target="_blank">outliers</a>, observations that have essentially different nature than most other observations. Their presence can lead to untrustworthy conclusions. The most complicated part of this task is <strong>to define a notion of “outlier”</strong>. After that, it is straightforward to identify them based on given data.</p>
<p>There are many techniques developed for outlier detection. Majority of them deal with numerical data. This post will describe the most basic ones with their application using <a rel="nofollow" href="http://dplyr.tidyverse.org/" target="_blank">dplyr</a> and <a rel="nofollow" href="https://echasnovski.github.io/ruler/" target="_blank">ruler</a> packages.</p>
<p>After reading this post you will know:</p>
<ul>
<li>Most basic outlier detection techniques.</li>
<li>A way to implement them using <code>dplyr</code> and <code>ruler</code>.</li>
<li>A way to combine their results in order to obtain a new outlier detection method.</li>
<li>A way to discover notion of “diamond quality” without prior knowledge of this topic (as a happy consequence of previous point).</li>
</ul>
</div>
<div>
<h1>Overview</h1>
<p>We will perform an analysis with the goal to find not typical diamonds listed in <a rel="nofollow" href="http://ggplot2.tidyverse.org/reference/diamonds.html" target="_blank">diamonds</a> dataset from <a rel="nofollow" href="http://ggplot2.tidyverse.org/" target="_blank">ggplot2</a> package. Here one observation represents one diamond and is stored as a row in data frame.</p>
<p>The way we will do that is <strong>by combining different outlier detection techniques to identify rows which are “strong outliers”</strong>, i.e. which might by considered outliers based on several methods.</p>
<p>Packages required for this analysis:</p>
<pre><code>library(dplyr)
library(tidyr)
library(ggplot2)
library(ruler)</code></pre>
</div>
<div>
<h1>Outlier detection methods</h1>
<p>To do convenient outlier detection with <code>ruler</code> it is better to define notion of <strong>non-outlier</strong> in form of the rule “Observation <strong>is not</strong> an outlier if …”. This way actual outliers are considered as <strong>rule breakers</strong>, objects of interest of <code>ruler</code> package. <strong>Note</strong> that definition of non-outlier is essentially a definition of outlier because of total two possibilities.</p>
<div>
<h2>Z-score</h2>
<p><a rel="nofollow" href="https://en.wikipedia.org/wiki/Standard_score" target="_blank">Z-score</a>, also called a standard score, of an observation is [broadly speaking] a distance from the population center measured in number of normalization units. The default choice for center is sample mean and for normalization unit is standard deviation.</p>
<p>⬛ <em>Observation is not an outlier based on z-score</em> if its absolute value of default z-score is lower then some threshold (popular choice is 3).</p>
<p>Here is the function for identifying non-outliers based on z-score:</p>
<pre><code>isnt_out_z &lt;- function(x, thres = 3, na.rm = TRUE) {
  abs(x - mean(x, na.rm = na.rm)) &lt;= thres * sd(x, na.rm = na.rm)
}</code></pre>
<p>It takes a numeric vector as input and returns logical vector of the same length indicating whether input value is a non-outlier.</p>
</div>
<div>
<h2>Z-score with MAD</h2>
<p><a rel="nofollow" href="https://en.wikipedia.org/wiki/Median_absolute_deviation" target="_blank">Median Absolute Deviation</a> is a robust normalization unit based on median as a population center. In order to use MAD “as a consistent estimator for the estimation of the standard deviation” one takes its value multiplied by a factor. This way base R function <a rel="nofollow" href="https://www.rdocumentation.org/packages/stats/versions/3.4.3/topics/mad?tap_a=5644-dce66f&tap_s=10907-287229" target="_blank">mad</a> is implemented.</p>
<p>⬛ <em>Observation is not an outlier based on MAD</em> if its absolute value of z-score with median as center and MAD as normalization unit is lower then some threshold (popular choice is 3).</p>
<pre><code>isnt_out_mad &lt;- function(x, thres = 3, na.rm = TRUE) {
  abs(x - median(x, na.rm = na.rm)) &lt;= thres * mad(x, na.rm = na.rm)
}</code></pre>
</div>
<div>
<h2>Tukey’s fences</h2>
<p><a rel="nofollow" href="https://en.wikipedia.org/wiki/Outlier#Tukey&amp;%2339;s_fences" target="_blank">Tukey’s fences</a> is a technique used in box plots. The non-outlier range is defined with <span>\([Q_1 – k(Q_3 – Q_1),~ Q_3 + k(Q_3 – Q_1)]\)</span>, where <span>\(Q_1\)</span> and <span>\(Q_3\)</span> are the lower and upper quartiles respectively, <span>\(k\)</span> – some nonnegative constant (popular choice is 1.5).</p>
<p>⬛ <em>Observation is not an outlier based on Tukey’s fences</em> if its value lies in non-outlier range.</p>
<pre><code>isnt_out_tukey &lt;- function(x, k = 1.5, na.rm = TRUE) {
  quar &lt;- quantile(x, probs = c(0.25, 0.75), na.rm = na.rm)
  iqr &lt;- diff(quar)
  
  (quar[1] - k * iqr &lt;= x) &amp; (x &lt;= quar[2] + k * iqr)
}</code></pre>
</div>
<div>
<h2>Mahalanobis distance</h2>
<p>All previous approaches were created for univariate numerical data. To detect outliers in multivariate case one can use <a rel="nofollow" href="https://en.wikipedia.org/wiki/Mahalanobis_distance" target="_blank">Mahalanobis distance</a> to reduce to univariate case and then apply known techniques.</p>
<p>⬛ <em>Observation is not an outlier based on Mahalanobis distance</em> if its distance is not an outlier.</p>
<pre><code>maha_dist &lt;- . %>% select_if(is.numeric) %>%
    mahalanobis(center = colMeans(.), cov = cov(.))

isnt_out_maha &lt;- function(tbl, isnt_out_f, ...) {
  tbl %>% maha_dist() %>% isnt_out_f(...)
}</code></pre>
<p>This function takes as input a data frame of interest (with possible non-numeric columns which are ignored) and function performing univariate outlier detection. It returns a logical vector of the same length as number of rows in input data frame.</p>
<p>To read more about practical usefulness of Mahalanobis distance in detecting outliers go to <a rel="nofollow" href="https://www.steffenruefer.com/2016/12/outlier-detection-with-mahalanobis-distance/" target="_blank">Steffen’s very helpful post</a>.</p>
</div>
</div>
<div>
<h1>Using dplyr and ruler</h1>
<div>
<h2>Definition of non-outlier row</h2>
<p>Package <code>ruler</code>, based on <code>dplyr</code> grammar of data manipulation, offers tools for validating the following data units: data as a whole, group [of rows] as a whole, column as a whole, row as a whole, cell. Our primary interest is row as a whole. However, using this framework, we can construct several approaches for definition of the non-outlier row:</p>
<ol>
<li><em>Row is not an outlier based on some column</em> if it doesn’t contain outlier (computed based on the target column) on the intersection with that column. In other words, first a univariate outlier detection is performed based solely on data from target column and then all rows containing non-outliers are named non-outlier rows.</li>
<li><em>Row is not an outlier based on Mahalanobis distance</em> if its distance (computed based on the selected numeric columns) is not an outlier.</li>
<li><em>Row is not an outlier based on grouping</em> if it is a part of a non-outlier group [of rows]. A group [of rows] is not an outlier if its summary value is not an outlier among summary values of other groups.</li>
</ol>
<p><strong>Note</strong> that all listed approached depend on the choice of the univariate outlier detection method. We will use all three previously listed univariate techniques.</p>
<pre><code>isnt_out_funs &lt;- funs(
  z = isnt_out_z,
  mad = isnt_out_mad,
  tukey = isnt_out_tukey
)</code></pre>
</div>
<div>
<h2>Implementation</h2>
<p>In <code>ruler</code> framework rules are defined in <a rel="nofollow" href="https://echasnovski.github.io/ruler/reference/rule-packs.html" target="_blank">packs</a> (to learn more go to <code>ruler</code> <a rel="nofollow" href="https://echasnovski.github.io/ruler/" target="_blank">README</a> and <a rel="nofollow" href="https://echasnovski.github.io/ruler/articles/" target="_blank">vignettes</a>).</p>
<div>
<h3>Column based non-outlier rows</h3>
<p>For <code>diamonds</code> dataset rules for column based non-outlier rows can be defined based on 7 numeric columns and 3 presented univariate detection methods. There is a convenient way of computing all them at once using <a rel="nofollow" href="http://dplyr.tidyverse.org/reference/summarise_all.html" target="_blank">scoped variant</a> of <code>dplyr::transmute()</code>:</p>
<pre><code>diamonds %>% transmute_if(is.numeric, isnt_out_funs)
## # A tibble: 53,940 x 21
##   carat_z depth_z table_z price_z   x_z   y_z   z_z carat_mad depth_mad
##     </code></pre></div></div></div>