# How to find the optimal cut-off point to minimize both the FNR and FPR in R?

I should find the optimal threshold to minimize both the false positive rate and false negative rate. An equal weight between these two rates should be assumed. I write the following code:

library(ROCR)
pred=prediction(data$$decile_score/10, data$$two_year_recid)
perf=performance(pred, measure="fnr",x.measure="fpr")

opt.cut = function(perf, pred)
{
cut.ind = mapply(FUN=function(x, y, p){
d = (x - 0)^2 + (y-1)^2
ind = which(d == min(d))
c(False_negative_rate = 1-y[[ind]], False_positive_rate = x[[ind]],
cutoff = p[[ind]])
}, [email protected], [email protected], pred@cutoffs)
}

print(opt.cut(perf, pred))

It throws out this result:

[,1]
False_negative_rate    0
False_positive_rate    0
cutoff               Inf

However, I think there is something wrong with my code.

• I like to use pacman, and p_load, because if it isn't installed it gets installed. Oct 26, 2022 at 13:30

Here are the libraries I used:

tidyverse,
reshape2)

url("https://raw.githubusercontent.com/propublica/compas-analysis/master/compas-scores-two-years.csv"),
sep=",")

Your data has 7214 rows and 53 columns, but you are only looking at 2. They are decile_score and two_year_recid.

Here is how I check the type, so I know they can plot well together.

data %>% select(decile_score, two_year_recid) %>% str()

This is the result:

'data.frame':   7214 obs. of  2 variables:
$$decile_score : int 1 3 4 8 1 1 6 4 1 3 ...$$ two_year_recid: int  0 1 1 0 0 0 1 0 0 1 ...

I looks like two numbers. Looks like the bottom is binomial, and the top is integers 1 to 10.

Here is a plot of the elements of prediction, to start to make sure that whatever "prediction" gives isn't junk.

ggplot(data, aes(x=decile_score/10, y=two_year_recid)) +
geom_point() +
geom_smooth(method="glm", formula = 'y ~ x ')

Here is the plot it yields:

Yes, nearly everything is wrong with the plot, but it still helps us ask questions.

In a very general sense, as decile score goes up, rate of recidivism goes up. As decile score goes down, rate of recidivism goes down. The place where recidivism chance is 50/50 looks to be aproximately the 50th percentile, or maybe a little above it.

Instead of using your pred, I'm going to use a proper linear model, and I'm going to turn recidivism into a proper data type. It is currently a float, but it should be a boolean, and the fit and prediction functions can respond to that.

df <- data.frame(x=data$$decile_score/10, y=data$$two_year_recid %>% as.factor)

est_glm <- glm(y~x, data=df, family = binomial())
summary(est_glm)

This yields:

The formula being fit in logistic regression is this:

The place of 50/50 odds occurs when:
$$\beta_0 + \beta_1 \cdot x = 0$$

solving that for $$x$$ gives:
$$-1.408958 + 2.650818 \cdot x = 0$$
$$2.650818 \cdot x = 1.408958$$
$$x = \frac {1.408958} {2.650818}$$

which evaluates to 53.15183%

This suggests that whatever decile is, it slightly under-predicts recidivism.

Some of the diagnostic plots from the glm fit can be interesting.

This suggests that the linear model in the exponent might not be the best. (I would test a cubic, and compare AIC values)

This one says there may be some outliers, a few points that corrupt the whole analysis.

Here is a violin plot that might speak to the nature of the population.

Why are there so many more non-recidivists at lower decile than high decile when the population of recidivists is nearly constant across deciles? There is a great reason for it, but without reading the paper, or understanding the sampling this can look like a surprising result.

So now for the tnr/tpr parts.

Here is the code I used:

cut_list <- seq(from=min(Preds), to=max(Preds), length.out=101)
store <- matrix(NA, nrow=length(cut_list), ncol=3) %>% as.data.frame()
names(store) <- c("cutoff", "tnr", "tpr")

for(i in 1:length(cut_list)){

preds2 <- ifelse(Preds > cut_list[i],1,0) %>% as.vector() %>% as.factor()

tp_num <- length(which(preds2 ==1 & df$$y == 1)) tp_den <- length(which(preds2 ==1 & df$$y == 1)) + length(which(preds2 ==0 & df$y == 1)) tpr <- tp_num/tp_den tn_num <- length(which(preds2 ==0 & df$$y == 0)) tn_den <- length(which(preds2 ==0 & df$$y == 0)) + length(which(preds2 ==1 & df$y == 0))

tnr <- tn_num/tn_den

store[i,1] <- cut_list[i]
store[i,2] <- tnr
store[i,3] <- tpr

}

I first plot each like this:

store2 <- melt(store, id.vars=1)

ggplot(store2, aes(x=cutoff, y=value, color=variable)) +
geom_point() + geom_path()

which gives this:

There are many ways to consider these values. Optimal means you have a measure of best, but I don't see one clearly shown here. If you consider them equal, you can get something like this.

Here is the code:

store3 <- store %>% mutate(both=tnr+tpr)

ggplot(store3, aes(x=cutoff, y=both)) +
geom_point() + geom_path() +
geom_vline(xintercept = 0.415, color="Red") +
geom_text(data = NULL, x = 0.415-0.025, y = 1.05, label = "x=0.415") +
geom_vline(xintercept = 0.478, color="Green") +
geom_text(data = NULL, x = 0.478+0.025, y = 1.05, label = "x=0.478") +
geom_vline(xintercept = 0.5*0.478+0.5*0.415, color="Black")

Here is the plot:

While you can get the same (big-chunk) result between 41.5% and 47.8%, you are likely to get the best result around 44.65%, which suggests that a likelihood above decile score above 44.65% is more likely than not to have recidivism within 2 years.

Like I said, I think the linear model has problems, and the data has a few points doing more of the work than they should. I don't trust someone else's model if I can't understand it, because that's a really fast way to make a big error, and I don't know where decile_score came from in detail.