Predict probability when model was trained in balanced dataset

Question

I have a dataset of about 1M observation and I had to predict a response that occurs only about 10.000 times (1%).

I decided to train a random forest, but this takes a lot of time to train because the data is too large for my hardware. So I decided to take a sample, but an aleatory sample would be already too large to have a minimum quantity of response. (If I take 10% aleatory, i would have only 1000 response)

Then I took a stratified sample. All responses and 10.000 aleatory non-responses, and trained my model in this dataset.

But now I need to rescale the probability so I have the real probability of the observation to be response.

I tried to simulate this problem with this code in R. Training a model in a balanced dataset and another one in the unbalanced data. But those models are not very correlated and I didn't find a good way to tranform the probabilities to the original unbalanced scale.

---

Found this is a good reference for future readers

---

I found that for logistic regression I can do this by just changing the intercept this way:

$$ \hat{\beta_0} = \hat{\beta_0^*} - log(\frac{\gamma_1}{\gamma_2})$$

Where $ \gamma_1 = Pr(Z=1|Y=1)$ and $ \gamma_2 = Pr(Z=1|Y=0)$. $Z$ is the an aleatory variable indicating if the observation is in the reduced dataset.

This was found in this book (in portuguese) page 216

---

simulate_data <- function(n){
  X <- data.frame(matrix(runif(n*20), ncol = 20))
  list(
    X = X,
    Y = rbinom(n, size = 1, prob = apply(X, 1, sum) %>% pnorm(mean = 13)) %>% as.factor()  
  )
}

balance <- function(X, Y){
  X <- rbind(
    X[Y == "1",],
    X[Y == "0",] %>% sample_n(length(Y[Y == "1"]))
    )
  return(list(
    X = X,
    Y = as.factor(c(rep(c(1,0), each = length(Y[Y == "1"]))))
  ))
}

train <- simulate_data(100000)

library(randomForest)
m_desb <- mean(train$Y == "1")
modelo_desb <- randomForest(train$X,train$Y, ntree = 200, cutoff = c(1-m_desb, m_desb), nodesize = 30, mtry = 8)
bal <- balance(train$X, train$Y)
m_bal <- mean(bal$Y == "1")
modelo_bal <- randomForest(bal$X,bal$Y, ntree = 100, cutoff = c(1-m_bal, m_bal), nodesize = 50)

---

Code for plot

library(ggplot2)
data.frame(
  unbalanced = predict(modelo_desb, newdata = test$X, type = "prob")[,2],
  balanced = predict(modelo_bal, newdata = test$X, type = "prob")[,2]
) %>%
  ggplot(aes(x = balanced, y = unbalanced)) +
  geom_point(size = 0.3) +
  xlim(0,1) +
  geom_smooth() +
  geom_hline(yintercept = m_desb, linetype = "dashed") +
  geom_vline(xintercept = m_bal, linetype = "dashed")

Answer 1

Daniel, What you did goes under the name of "oversampling." There is a sample of some "real" population, and you replace it with a sample from a "manufactured" population. The problem that makes sense in application is the estimation of

$$P_r(Y=1|X) = \text{probability of response=1 in the $\mathbf{real}$ population given the predictor $X$} $$

but by using an oversample you are estimating $$P_m(Y=1|X) = \text{probability of response=1 in the $\mathbf{manufactured}$ population given the predictor $X$}$$

The two probabilities are related. I'll worked the details. I'll pretend the predicor $X$ is discrete. If $X$ takes numerical values one has to replace some probabilities by probability densities. $$\dots\dots\dots$$

To simplify the notation, let $\pi_1 = P_r(Y=1)$ and $\mu_1 = P_m(Y=1)$ be the probabilities of response in the real and manufactured populations, let $$ L_r = \frac{P_r(X=x|Y=1)}{P_r(X=x|Y=0)} = \frac{\frac{P_r(Y=1|X=x)}{P_r(Y=0|X=x)}}{\frac{\pi_1}{1-\pi_1} }$$ be the odds ratio of $Y=1$, i.e.: the ratio of the odds among cases with $X=x$ and the odds in the general $\mathbf{real}$ population. Finally, let $L_m$ be the corresponding ratio in the $\mathbf{manufactured}$ population.

By Bayes' Theorem: $$ P_r(Y=1|X=x) = \frac{P_r(Y=1,X=x)}{P_r(X=x)} = \\ =\frac{P_r(X=x|Y=1)\space \pi_1}{P_r(X=x|Y=1)\space\pi_1 + P_r(X=x|Y=0)\space (1 - \pi_1)} = \\ =\frac{L_r\space \pi_1}{L_r\space\pi_1 + \space (1 - \pi_1)} \tag{1} $$ In a similar way, we get an analogous result for the manufactured population: $$ P_m(Y=1|X=x) = \frac{L_m\space \mu_1}{L_m\space\mu_1 + \space (1 - \mu_1)} \tag{2} $$ Since the manufactured sample is a random sample, stratified by $Y$, the conditional distribution of X within responders is the same as in the real population. Same as for non responders, i.e.:
$$ P_r(X=x|Y=j) = P_r(X=x|Y=j) $$ for $j=0,1.$ If the sample stratified by values of Y was anything other than random sample these would not be true. It follows that $\boxed{ L_r = L_m }$. Next we solve for $L_m$ in terms of $P_m(Y=1|X)$ from (2) and replace in (1).

$$\dots\dots\dots$$ $\mathbf{Digression}$: Here is my easy way to carry the steps, without mess. Two non-zero vectors $\mathbf{v_1}$, $\mathbf{v_2}$ are parallel iff there is $\lambda \ne 0$ such that $\mathbf{v_1}=\lambda \mathbf{v_2}.$ Below I will use this idea, and I will not care about the exact value of $\lambda$, so I will be using "$\lambda$" as a short-hand for "$\mathbf{\text{some non-zero mumber}}$." $\mathbf{\text{End Digression}}$: $$\dots\dots\dots$$

The easy way to solve is to observed that for non-zero messy values of $\lambda$ (not the same in each occurrence!) one has:

$$ \begin{bmatrix} P_r(Y=1|X) \\ 1 \\ \end{bmatrix} = \lambda \begin{bmatrix} \pi_1 &0\\ \pi_1 &1-\pi_1 \end{bmatrix} \begin{bmatrix} L_r\\ 1 \end{bmatrix} , $$ and $$ \begin{bmatrix} P_m(Y=1|X) \\ 1 \\ \end{bmatrix} = \lambda \begin{bmatrix} \mu_1 &0 \\ \mu_1 &1-\mu_1 \end{bmatrix} \begin{bmatrix} L_m\\ 1 \end{bmatrix} . $$

Therefore,

$$ \begin{bmatrix} L_m\\ 1 \end{bmatrix} = \lambda \begin{bmatrix} \mu_1 &0 \\ \mu_1 &1-\mu_1 \end{bmatrix}^{-1} \begin{bmatrix} P_m(Y=1|X) \\ 1 \\ \end{bmatrix} , $$

and so (remember that here $\lambda$ stands for "some non-zero number")

$ \space \begin{bmatrix} P_r \\ 1 \end{bmatrix} = \lambda \begin{bmatrix} \pi_1 &0 \\ \pi_1 &1-\pi_1 \end{bmatrix} \begin{bmatrix} L_r \\ 1 \end{bmatrix} = \\ \text{ }=\lambda \begin{bmatrix} \pi_1 &0 \\ \pi_1 &1-\pi_1 \end{bmatrix} \begin{bmatrix} \mu_1 &0 \\ \mu_1 &1-\mu_1 \end{bmatrix}^{-1} \begin{bmatrix} P_m \\ 1 \\ \end{bmatrix} = \\ \text{ }=\lambda \begin{bmatrix} \pi_1 (1- \mu_1) &0 \\ \pi_1 - \mu_1 & \mu_1 (1- \pi_1) \end{bmatrix} \begin{bmatrix} P_m \\ 1 \\ \end{bmatrix} = \lambda \begin{bmatrix} \pi_1 (1-\mu_1) P_m \\ (\pi_1 - \mu_1) \; P_m + \mu_1 (1- \pi_1) \end{bmatrix}. $
Thus, $$ P_r = \frac{\pi_1 (1- \mu_1) P_m}{(\pi_1 - \mu_1) \; P_m + \mu_1 (1- \pi_1) } $$

$$\dots\dots\dots$$ Example: Let's work the details of a Binomial model, $$P_m(Y=1|X) = \frac{e^{\beta_0 + \beta X}}{1+e^{\beta_0 + \beta X}} $$ or in the "where $\lambda$ is some non-zero scalar" notation (I would not had digressed before if I did not had ulterior motive.. :) ):

$$ \begin{bmatrix} P_m \\ 1 \\ \end{bmatrix} = \lambda \begin{bmatrix} e^{\beta_0 + \beta X} \\ 1 + e^{\beta_0 + \beta X} \end{bmatrix} $$

What is the implied model in the real population?

$ \space \begin{bmatrix} P_r \\ 1 \end{bmatrix} = \lambda \begin{bmatrix} \pi_1 (1- \mu_1) &0 \\ \pi_1 - \mu_1 &\mu_1 (1- \pi_1) \end{bmatrix} \begin{bmatrix} P_m \\ 1 \end{bmatrix} = $

$\;$ $ = \lambda \begin{bmatrix} \pi_1 (1- \mu_1) &0 \\ \pi_1 - \mu_1 &\mu_1 (1- \pi_1) \end{bmatrix} \begin{bmatrix} e^{\beta_0 + \beta X} \\ 1 + e^{\beta_0 + \beta X} \end{bmatrix}= $

$\;$ $ = \lambda \begin{bmatrix} \pi_1 (1- \mu_1) e^{\beta_0 + \beta X} \\ \pi_1 (1- \mu_1) e^{\beta_0 + \beta X} + \mu_1 (1- \pi_1) \end{bmatrix} = \lambda \begin{bmatrix} \frac{\pi_1 (1- \mu_1)}{\mu_1 (1- \pi_1)} e^{\beta_0 + \beta X} \\ 1 + \frac{\pi_1 (1- \mu_1)}{\mu_1 (1- \pi_1)} e^{\beta_0 + \beta X} \end{bmatrix} .$

If we let $\tau = \ln(\frac{\pi_1 (1- \mu_1)}{\mu_1 (1- \pi_1)})$, we can absorve this constant in the exponents to get: $$ \begin{bmatrix} P_r \\ 1 \end{bmatrix} = \lambda \begin{bmatrix} e^{\tau + \beta_0 + \beta X} \\ 1 + e^{\tau + \beta_0 + \beta X} \end{bmatrix} .$$
Taking the ratio and simplifying the non-zero constant in numerator and denominator we get that fitting a logistic model to the manufactured population results in an implied logistic model for the real population, $\mathbf{\text{with the same coefficients for X}}$ and with a difference in the constant (in the logistic model) given by: $$ \beta_{real} = \tau + \beta_0 $$

$$\dots$$ Note that, according to your reference, the ratio of $\gamma_1 = Pr(Z=1|Y=1)$ and $\gamma_0 = Pr(Z=1|Y=0)$ should come up. Indeed:
$$ \gamma_1 = Pr(Z=1|Y=1) = \frac{P(Z=1,Y=1)}{P(Y=1)} = \frac{P_r(Y=1|Z=1)P_r(Z=1)}{P_r(Y=1)} = \frac{P_m(Y=1)}{P_r(Y=1)} P_r(Z=1)= \frac{\mu_1}{\pi_1}P_r(Z=1) $$ likewise (i.e. change Y to 1-Y), $$ \gamma_0 = \frac{1-\mu_1}{1-\pi_1}P_r(Z=1) $$ so $$ \ln(\frac{\gamma_1}{\gamma_0}) = - \ln(\frac{\pi_1 (1-\mu_1)}{\mu_1 (1-\pi_1}) = \tau $$

$$\dots\dots\dots$$ Notes for full disclosure: I worked with the probability model. When one works with finite samples the example above suggests two ways of estimating the coefficients: * estimate coefficients using the sample from the real population * estimate coefficients using the manufactored populations

It terns out that this two estimators are not the same (it is obvious if one consideres one estimator is based on more cases than the other). Both estimators are asymtopically consistent, but it can be shown the one based on the manufactored population is more biased (forgot the reference :( ).

In the data science space we are more concern with the quality of the predictions than the parameters of the parameters used to make those predictions, so as long as you check results properly (e.g.: using a testing set to build models and another to validate them), the bias in the parameters should not deter us from using oversampling.

$$\dots\dots\dots$$

Answer 2

Try using Platt's trick for probability learning. This is just a logistic regression of the model output to the real probability (see paper )

Answer 3

I believe for your case a better optimization metric like F1 score will work better since you have unbalanced labels.