# Question

What is the real-life example of the benefit and application of the benefit of Bayesian regression?

Having read the items and it looks having the range of inference (possible values and likelihood) available is the benefit.

But how it will be utilised in real life (finance, engineering, ...)? I suppose without real-life applications of the benefit, it would not be so useful.

It seems both linear regression and Bayesian regression can produce similar predictions as below.

According to 3, the predictive distribution can give the confidence on the prediction if it is within the dense-color area because of the data is dense, but not in sparse area, eg. prediction at x=5 may not be trustworthy. But I suppose there should be more than that.

What does "having statistical inference" make a difference and how it is utilised in real life? Or Is it a matter of choice to use linear or Bayesian?

# Code

Taken from (3) and made a few changes. I am not completely familiar with the logic so if there are mistakes, please correct.

import math
import numpy as np
import matplotlib.pyplot as plt

# --------------------------------------------------------------------------------
# Gaussian basis function
# --------------------------------------------------------------------------------
def gaussian(mean, sigma):
"""
Args:
mean:
sigma:
"""
def _gaussian(x):
return np.exp(-0.5 * ((x - mean) / sigma) ** 2)
return _gaussian

# --------------------------------------------------------------------------------
# Design matrix
# --------------------------------------------------------------------------------
def phi(f, x):
bias = np.array([1])  # bias parameter(i = 0)
# bias+basis
return np.append(bias, f(x))

# --------------------------------------------------------------------------------
# Data generation utility
# --------------------------------------------------------------------------------
from numpy.random import rand
def uniform_variable_generator(samples):
_random = rand(samples)
return _random

def noise_generator(samples, mu=0.0, beta=0.1):
noise = np.random.normal(mu, beta, samples)
return noise

def sigmoid(x):
return 1 / (1 + np.exp(-x))

def generator_t_function(x):
#return np.sin(x)
return sigmoid(x)

def generator_X_function(x):
return 2 * np.pi * x
#return 2 * np.pi * x

# --------------------------------------------------------------------------------
# Observations
# --------------------------------------------------------------------------------
#X = np.array([0.02, 0.12, 0.19, 0.27, 0.42, 0.51, 0.64, 0.84, 0.88, 0.99])
#t = np.array([0.05, 0.87, 0.94, 0.92, 0.54, -0.11, -0.78, -0.79, -0.89, -0.04])
samples = 20

#X = np.array([0.02, 0.12, 0.19, 0.27, 0.42, 0.51, 0.64, 0.84, 0.88, 0.99])
#t = np.array([0.05, 0.87, 0.94, 0.92, 0.54, -0.11, -0.78, -0.79, -0.89, -0.04])
X = generator_X_function(uniform_variable_generator(samples))
t = generator_t_function(X) + noise_generator(samples, beta=0.1)

MAX_X = max(X)
NUM_X = len(X)
MAX_T = max(t)
NUM_T = len(t)

# --------------------------------------------------------------------------------
# Gaussian basis function parameters
# --------------------------------------------------------------------------------
sigma = 0.1 * MAX_X

# mean of gaussian basis function (11 dimension w1, w2, ... w11)
mean = np.arange(0, MAX_X + sigma, sigma)

# Basis function
f = gaussian(mean, sigma)

# --------------------------------------------------------------------------------
# Design matrix
# --------------------------------------------------------------------------------
PHI = np.array([phi(f, x) for x in X])

#alpha = 0.1
#beta = 9.0
alpha = 0.5  # larger alpha gives smaller w preventing overfitting (0 -> same with linear regression)
beta = 5   # Small beta allows more variance (deviation)

Sigma_N = np.linalg.inv(alpha * np.identity(PHI.shape[1]) + beta * np.dot(PHI.T, PHI))

mean_N = beta * np.dot(Sigma_N, np.dot(PHI.T, t))

# --------------------------------------------------------------------------------
# Bayesian regression
# --------------------------------------------------------------------------------
xlist = np.arange(0, MAX_X, 0.01)
plt.title("Bayesian regression")
plt.plot(xlist, [np.dot(mean_N, phi(f, x)) for x in xlist], 'b')
plt.plot(X, t, 'o', color='r')
plt.show()

# --------------------------------------------------------------------------------
# Linear regression
# --------------------------------------------------------------------------------
# w for linear regression parameter
#w = np.linalg.solve(np.dot(PHI.T, PHI), np.dot(PHI.T, t))
# --------------------------------------------------------------------------------
l = 0.05
regularization = np.identity(PHI.shape[1])
w = np.linalg.solve(
np.dot(PHI.T, PHI) + (l * regularization),
np.dot(PHI.T, t)
)

xlist = np.arange(0, MAX_X, 0.01)
plt.title("Linear regression")
plt.plot(xlist, [np.dot(w, phi(f, x)) for x in xlist], 'g')
plt.plot(X, t, 'o', color='r')
plt.show()

# --------------------------------------------------------------------------------
# Predictive Distribution
# --------------------------------------------------------------------------------
def normal_dist_pdf(x, mean, var):
return np.exp(-(x-mean) ** 2 / (2 * var)) / np.sqrt(2 * np.pi * var)

return np.dot(x, np.dot(A, x))

xlist = np.arange(0, MAX_X, 0.01)
#tlist = np.arange(-1.5 * MAX_T, 1.5 * MAX_T, 0.01)
tlist = np.arange(
np.mean(t) - (np.max(t)-np.min(t)),
np.mean(t) + (np.max(t)-np.min(t)),
0.01
)
z = np.array([
normal_dist_pdf(tlist, np.dot(mean_N, phi(f, x)),
1 / beta + quad_form(Sigma_N, phi(f, x))) for x in xlist
]).T

plt.contourf(xlist, tlist, z, 5, cmap=plt.cm.coolwarm)
plt.title("Predictive distribution")
plt.plot(xlist, [np.dot(mean_N, phi(f, x)) for x in xlist], 'r')
plt.plot(X, t, 'go')
plt.show()

• Starting from the first image in your question you have lost me... How do the images and the attached python code relate to the question about general use? – Frans Rodenburg Apr 29 '19 at 1:18
• Related/duplicate: stats.stackexchange.com/q/8347/176202 – Frans Rodenburg Apr 29 '19 at 1:22
• @FransRodenburg, thanks for the follow up. Intended show Linear and Bayesian could do the same as in the image as far as I experimented as in the attached code. Which make me think there may not be much benefit if they could do the similar, and only difference I saw was the density plot. I suppose this density would make the difference. So I wonder this density is being utilised and if so how. – mon Apr 29 '19 at 1:23
• A related question: stats.stackexchange.com/q/252577/103153 – Lerner Zhang Apr 29 '19 at 2:24
• @FransRodenburg, kindly provide a case where you have used Bayesian regression in real life and the reason why Bayesian is good fit for the case? – mon May 1 '19 at 1:23

We use a Bayesian model in forecasting retail sales with the SAP Unified Demand Forecast (UDF). The Bayesian approach offers two advantages:

1. We can use priors. For instance, suppose your supermarkets stocks a new item, and we would like to forecast its first Christmas sales. We can simply use the average effect of similar items as a prior for the new item's Christmas effect. Which has the added benefit that the prior gets automatically updated once we have seen the new item's first Christmas sales, so both the prior and the item's own data influence the prediction for the second Christmas.

2. The priors regularize. Our model is heavily over-parameterized, with seasonality, day of week, trend, holidays and tons of promotion predictors, so regularization is hugely important to keep the predictions under control. (Note that the pictures you show in your question do not show such an overparameterized model.)

Aspect 1 is more important to our marketers and to business users, and aspect 2 is more important to me and to statistician or data scientist users.

Yes, there are ways to address both issues without a Bayesian model, e.g., by pooling for aspect 1, and the Lasso or Elastic Net for aspect 2. We just chose the Bayesian approach, and I find the updating particularly elegant.

So, yes, Bayesian regression is used in real life, and it has paid my salary for the last couple of years, along with the salaries of a number of my colleagues. And while making sure your supermarket does not run out of your favorite shampoo.