Means , medians, and modes provide measures of central tendency. That is, they are single numbers that sit in the center of a collection of values. But how is this "center" defined? This post will provide some intuition into how they are related, and how they are intertwined with loss functions in machine learning.
using Plots, LossFunctions, Distributions
plotlyjs();
srand(0xdead);
Perhaps the most common measure of central tendency is the arithmetic mean, also commonly referred to as the average or simply the mean. The mean is defined as the sum of observations divided by the number of observations, and it is quite closely related to the concept of squared error.
#Generate some data from a normal distribution with mean 3 and standard deviation 1
dist = Normal(3,1);
data = rand(dist, 100);
This data has a "bell-curve" shape that we all know and love. Visual inspection shows that the center of the distribution should be somewhere around 3
histogram(data)
#Given a centrality estimator (i.e. the mean), sums the error between that estimator and each data point
#as measured by the loss function.
function average_loss(est::Number, data::Vector{<: Number}, loss::Loss)
n = length(data)
#Crudely broadcast the estimator to act like an array
ests = [est for i in 1:n]
value(loss, data, ests, AvgMode.Mean())
end
function plot_average_loss(data::Vector{<:Number}, loss::Loss; minrange=minimum(data), maxrange=maximum(data))
avgerrorfunc = est -> average_loss(est, data, loss)
#Central tendency estimator should be between minimum and maximum of data
p = plot(avgerrorfunc, linspace(minrange, maxrange, 200), xlabel="Estimator Value", ylabel="Average Error")
end
#Also make functions to plot losses
function plot_loss(loss::DistanceLoss; minrange=-2, maxrange=2)
lossfunc = dist -> value(loss, dist)
plot(lossfunc, linspace(minrange, maxrange, 200), xlabel="Distance from true value", ylabel="Loss value")
end
function plot_loss(loss::MarginLoss; minrange=-2, maxrange=2)
lossfunc = agreement -> value(loss, agreement)
plot(lossfunc, linspace(minrange, maxrange, 200), xlabel="Agreement with true value", ylabel="Loss value")
end
First, let's plot the squared error loss function
plot_loss(L2DistLoss())
Similarly to the loss function itself, the average loss function over the dataset looks like a quadratic function.
plot_average_loss(data, L2DistLoss())
Here we can see that the average error as we change the "center" looks like a parabola centered around 3. We defined the normal distribution with mean 3, so let's see what the actual mean is
mean(data)
This is close to 3, but is it the actual minimizer? Let's take a closer look at values around the mean to confirm that the mean of the data minimizes the average squared error to every other point in the data.
#L2DistLoss is just another way to say "Squared Error"
plot_average_loss(data, L2DistLoss(), minrange=mean(data)-0.001, maxrange=mean(data)+0.001)
The true minimum of the average loss is in fact, the mean, as shown in this plot of the average loss around a small neighborhood.
Of course, not all data comes from a perfect normal distribution. Here, we'll see where other loss functions, and hence measures of central tendency can be valuable. Outlying values affect different centrality measures differently, and I'll try to show here that this is because they assume different error metrics.
data_outliers = rand(dist, 100);
#The first 10 entries are outliers
data_outliers[1:10] = rand(Normal(10, 5), 10);
histogram(data_outliers)
Now the data still looks kind of like a normal distribution, but there are some clear outliers. The eye can still identify that a kind of "center" is at 3, but the mean may tell something different.
function plot_average_loss(data::Vector{<:Number}, losses::Array; minrange=minimum(data), maxrange=maximum(data))
errfuncs = [(est -> average_loss(est, data, loss)) for loss in losses]
labels = [string(loss) for loss in losses]
#Central tendency estimator should be between minimum and maximum of data
p = plot(errfuncs, linspace(minrange, maxrange, 200),
xlabel="Estimator Value", ylabel="Average Error", label=labels)
end
function plot_losses(losses::Array; minrange=-2, maxrange=2)
lossfuncs = [dist -> value(loss, dist) for loss in losses]
plot(lossfuncs, linspace(minrange, maxrange, 200), xlabel="Agreement", ylabel="Loss value")
end
The absolute value loss function increases linearly on both sides.
plot_loss(L1DistLoss())
We can plot the average error for the squared error function (blue) and the absolute value function (red). It's pretty easy to see that the functions are minimized in different places. The value that minimizes the average absolute error is smaller than that for squared error. In fact, these are exactly the median and mean, respectively
A couple extra things to notice. First, the squared error is larger everywhere than the absolute error, and grows faster in both directions. This is one explanation as to why the mean is more sensitive to outliers than the median.
p = plot_average_loss(data_outliers, [L2DistLoss(), L1DistLoss()], minrange=0, maxrange=10)
vline!([mean(data_outliers)], label="mean")
vline!([median(data_outliers)], label="median")
Let's take a closer look at the error curves. First, we can see that the average L2 (squared error) loss is minimized at the mean, almost 4.0 .
plot_average_loss(data_outliers, L2DistLoss(),
minrange=mean(data_outliers)-0.1, maxrange=mean(data_outliers)+0.1)
mean(data_outliers)
However, the L1 (absolute-value) loss is minimized at the median, much closer to 3.
plot_average_loss(data_outliers, L1DistLoss(),
minrange=median(data_outliers)-0.1, maxrange=median(data_outliers)+0.1)
median(data_outliers)
Of course, more loss functions exist for regression, such as the L1 Epsilon-Insensitive loss (in Support Vector Regression), the Huber Loss function, and the log-cosh loss function (infinitely-differentiable Huber Loss). While these all share the robustness property with the absolute-value loss function, their average values aren't minimized at exactly the median. Let's take a look... (Absolute Value in blue, Huber in red, Epsilon Insensitive in green, Log-Cosh in purple)
plot_average_loss(data_outliers, [L1DistLoss(), HuberLoss(), L1EpsilonInsLoss(0.5), LogitDistLoss()],
minrange=2, maxrange=5)
vline!([median(data_outliers)], label="median")
Statistics other than the mean are used frequently in reporting things like income, where there are very few extremely high incomes that bias the mean upwards. Hence, household income is frequently reported in terms of median income.
These estimators also represent some of the simplest machine learning models. The mean is equivalent to a least-squares regression model with just an intercept, while the median is an L1-regression intercept. These are common baselines to use when fitting models to data, to make sure that model predictions are better than a constant.
So far, we've looked at the mean and median as estimates of central tendency in data. These are part of a class of estimators called M-estimators, developed by Huber in 1964 and 1967. M-estimators are scalar values that minimize the sum of loss functions between themselves and a sample of data. However, so far the data has been sampled from a continuous distribution, with loss functions that measure the difference between an estimate and a true value. These loss functions are normally used for regression problems, where the data can take on any value in a continuous range (i.e. temperature, blood pressure, power, etc.).
Often, data is more restricted in the values it can take. For example, it may be binary-valued, ordinal (discrete with distinct levels), or categorical (discrete without any ordering). For the binary case, common loss functions include the Logistic loss (in logistic regression) and the Hinge loss (in Support Vector Machines). The values that minimize these loss functions are less well-known, but constitute the bias term in logistic regression and SVMs, respectively.
#Generate data from a bernoulli distribution (i.e. an unbalanced coin flip)
#Failure is represented by -1, success is 1
srand(123)
data_bin = 2 .* rand(Bernoulli(0.65), 100) .- 1;
plot_loss(LogitMarginLoss())
Unlike the distance loss cases, the plot of average loss value with an estimator does not mirror a visualization of the loss function itself. Here I'll plot the logistic loss function. It decreases monotonically because it treats the x axis as a measure of confidence. A high confidence in a correct answer (high positive agreement) has low error, while high confidence in a wrong answer (negative agreement) has higher error.
plot_average_loss(data_bin, LogitMarginLoss(), minrange=0, maxrange=1.1)
vline!([mean(data_bin)], label="mean")
vline!([median(data_bin)], label="median")
mean(data_bin)
median(data_bin)
The M-estimator of the logistic loss function is different from both the mean and the median. While the mean is around 0.24, and the median is 1 (due to the class distribution), the minimum of the average loss is achieved in between. This is equivalent to the intercept in a logistic regression model with no predictors.
There are also other loss functions for classification, which generally decrease as agreement increases. Hinge loss and Squared Hinge loss (used in Support Vector Machines) are a couple of them. This plot shows the logistic loss (blue) alongside the hinge loss (red) and squared hinge loss (green).
plot_losses([LogitMarginLoss(), HingeLoss(), L2HingeLoss()])
Let's also look at the error achieved by an estimator for these loss functions (logistic loss in blue, hinge loss in red, and squared hinge loss in green)
plot_average_loss(data_bin, [LogitMarginLoss(), HingeLoss(), L2HingeLoss()],
minrange=0, maxrange=1.1)
vline!([mean(data_bin)], label="mean")
vline!([median(data_bin)], label="median")
Here, we see that different classification losses are minimized at different values. The hinge loss has its m-estimator at the median, while the squared hinge loss has its m-estimator at the mean. As above, the logistic loss is minimized for a constant in between the two.