First Example

Bag with unknown number of blue and white balls inside. Two marbles are drawn with replacement. Four possibilities are: ww, bb, wb, and bw.

Supoose we know that the expected no of blue marbles over two draws is exactly 1. This is the expected value constraint we apply.

We consider four prob dists.

dist <- list(
  A = c(1/4, 1/4, 1/4, 1/4),    # binomial with n = 2, p = 0.5
  B = c(2/6, 1/6, 1/6, 2/6),    
  C = c(1/6, 2/6, 2/6, 1/6), 
  D = c(1/8, 4/8, 2/8, 1/8)
)

withr::with_par(
  list(mfrow = c(2, 2)), 
  code = {
    for(i in seq_along(dist)){
      plot(1:4, dist[[i]], xaxt = "n", yaxt = "n", ylim = c(0, 0.6), 
           type ="b", axes = FALSE, xlab = "", ylab = "",
           lwd = 2, col = 'steelblue', pch = 16, )
      axis(side = 1, at = 1:4, labels = c("ww", "bw", "wb", "bb"), )
      text(1, 0.45, labels = names(dist)[i], cex = 2)
    }
  }
)

A is the binomial distribution with n = 2, and p = 0.5 (with bw and wb collapsed into same outcome type).

B, C and D are not binomial but have same expected value.

sapply(dist, function(p) sum(p * c(0, 1, 1, 2)))

## A B C D 
## 1 1 1 1

Entropy of each:

round(sapply(dist, entropy), 3)

##     A     B     C     D 
## 1.386 1.330 1.330 1.213

A has highest entropy. This is since H increases as probability distribution becomes more even.

In this example the expected value is 1, which is why the distribution over outcomes can be flat and remain consistent with the constraint.

Exponential Family

Every member of this family a max entropy distribution for some set of constraints.

Exponential distribution is constrained to be zero or positive. It is a dist of distance/duration, measurements that represent displacement from some point of reference (in time or space). If the prob of an event is constant in time or across space, then the dist of events tends towards exponential. It has the max entropy among all non-negative continuous distributions with the same average displacement.

Its shape is described by single parameter \(\lambda\), or the avg displacement \(\lambda^{-1}\).

Gamma Distribution - also constrained to be >= 0. Unlike exponential it can have a peak above zero. If an event can only happen after 2 or more exp distributed events happen, the resulting waiting times will be gamma distributed.

Has highest entropy among distributions with the same mean and same average logarithm. Shape is described by 2 params.

Poisson Distribution - a count distribution like binomial. Mathematically a special case of binomial, when \(n\) is very large and \(p\) is very small. Parameter \(\lambda=np\). It is used for counts that never get close to any theoretical maximum. Described only by \(\lambda\), the rate of events.

Linking linear models to distributions

To build a regression model from any of the exp family distributions involves creating linear model(s) and attaching them to the parameter(s) of the distribution that describe its shape. A link function is often required for this. A link function’s job is to map the linear space of a model onto the non-linear space of the paramater. Two common links are logit link and log link.

Logit Link is for mapping a parameter that is defined as a probability mass (thus, must lie between 0 and 1).

$$ \[\begin{align} \operatorname{logit}(p_i) &= \log \frac{p_i}{1 - p_i} \tag{Log-Odds} \\ \log \frac{p_i}{1 - p_i} &=\alpha + \beta x_i \\ \text{Arranging terms we get} \implies \\ p_i &= \frac{\exp (\alpha + \beta x_i)}{1 + \exp (\alpha + \beta x_i)} \end{align}\]

$$

Above function is called logistic. Also the Inverse Logit, since it inverts the logit transform. When we use a logit transform, we are defining the distribution’s parameter value to be the logistic transform of the linear model.

x <- seq(-3, 3, length.out = 100)
log_odds <- 2 * x
p <- exp(log_odds) / (1 + exp(log_odds))

withr::with_par(
  list(mfrow = c(1, 2)), 
  code = {
    plot(x, log_odds, type = "l", ylab = "log-odds", col = 'navyblue')
    
    plot(x, p, type = "l", ylab = "probability", col = 'navyblue')
    abline(h = c(0, 1), lty = 2)
    
  }
)

This affects interpretation, since a unit change in predictor does not correspond with a constant change in the mean of the outcome variable. (Comparing to interaction terms, every parameter sort of interacts with itself now, becuase the impact of a change in a prector depends upon the value of the predictor before change).

Second common link is LOG LINK, which maps a parameter defined only over positive real values onto a linear model. For instance sd of a Gaussian distribution can only be positive (neither negative or zero).

$$ \[\begin{align} y_i &\sim Normal(\mu, \sigma_i) \\ \log (\sigma_i) &= \alpha + \beta x_i \\ \text{This implies that} \\ \sigma_i &= exp(\alpha + \beta x_i) \end{align}\] $$ Thus, log link assumes parameter’s value is the exp of the linear model, and impolies exponential scaling of the outcome with the predictor variable.

x <- seq(-2, 2, length.out = 100)
y <- 2 * x
ly <- exp(y)

withr::with_par(
  list(mfrow = c(1, 2)), 
  code = {
    plot(x, y, type = "l", ylab = "log measurement", col = 'navyblue')
    
    plot(x, ly, type = "l", ylab = "original measurement", col = 'navyblue')
    
  }
)