# Predict probability when model was trained in balanced dataset

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")

• Just to clarify, you want to find probabilities of a response that only occurs 1% of the time using logistic regression? Or do you still want to use random forest? Also, did the 10% sample do poorly in cross-validation? Jul 28, 2016 at 15:18
• @Hobbes I want to keep using random forests. Yes, the classification acuracy in test set using the balanced dataset was 20% better. Jul 28, 2016 at 15:21
• Have you considered a different probability threshold for classification purposes (maybe everything p>.10 instead of p>.50)? I would recommend also examining precision and recall in a situation like this (including ROC/AUC/F1/Brier Score). In my experience trees don't bode to well with major imbalance and trying to turn a classification result into a modeled probability is a tough task. Jul 28, 2016 at 18:29
• I used the ROC curve to choose the probability threshold for both cases in my simulations. And I examined the KS. Jul 28, 2016 at 18:34

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$$

• Thank you very much for your answer! I learned a lot with it. But I still can't solve my problem. I understand the binomial model I can find new estimates for the coefficients using this approach. But when I'm using random forest is it still possible practically? Aug 2, 2016 at 17:53
• Sure. The computation is about over-sampling, and how to relate the model in the manufacured sample and in the real sample. If works for binary variables, I've never use it in other case, so I have not thought of the corresponding correction in cas of multinomial response (trees are great to model multinomial targets). Aug 5, 2016 at 0:36

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

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