Agathe Fernandes Machado
email: fernandes_machado.agathe@courrier.uqam.ca

Arthur Charpentier

Emmanuel Flachaire

Ewen Gallic

François Hu

From Uncertainty to Precision: Enhancing Binary Classifier Performance through Calibration

Keywords: Calibration, Binary classification, Local regression

The assessment of binary classifier performance traditionally centers on discriminative ability using metrics, such as accuracy. However, these metrics often disregard the model's inherent uncertainty, especially when dealing with sensitive decision-making domains, such as finance or healthcare. Given that model-predicted scores are commonly seen as event probabilities, calibration is crucial for accurate interpretation. In our study, we analyze the sensitivity of various calibration measures to score distortions and introduce a refined metric, the Local Calibration Score. Comparing recalibration methods, we advocate for local regressions, emphasizing their dual role as effective recalibration tools and facilitators of smoother visualizations. We apply these findings in a real-world scenario using Random Forest classifier and regressor to predict credit default while simultaneously measuring calibration during performance optimization.

Introduction

Binary classification tasks are prevalent in learning algorithms, as diverse scenarios require binary decisions. Examples include predicting default risk or accident occurrence in insurance or finance as well as disease likelihood in healthcare. To improve reliability, particularly in sensitive decision-making contexts, a classifier must possess strong discriminatory capabilities. Typically, classifiers are trained to optimize goodness-of-fit criteria, often based on the accuracy of class predictions. However, goodness-of-fit criteria, such as accuracy or AUC, do not consider the varying confidence levels assigned by the algorithm to each prediction. If the sole objective is effective class prediction, then the classifier fulfills its purpose. Nevertheless, there are instances where interest extends beyond the predicted class to the associated likelihood. This occurs when predicting loan repayment defaults [START_REF] Liu | Applying machine learning algorithms to predict default probability in the online credit market: Evidence from china[END_REF] or accident incidences, as risk transfer pricing is usually tied directly to event probabilities. In such cases, the model-predicted scores of classifiers are often interpreted as event probabilities. Yet, in these examples, achieving accurate interpretation as probabilities necessitates effective calibration of the model. A well-calibrated model ensures precise understanding of the predicted scores as probabilities. For instance, if a model assigns a predicted probability of 80% to events, the observed proportion of those events occurring over the long runaccording to [START_REF] Dawid | The well-calibrated bayesian[END_REF] terminology-should ideally align with the predicted value of 80%.

Simple classifiers such as Logistic Regression models typically exhibit overall calibration [START_REF] Mildenhall | A systematic relationship between minimum bias and generalized linear models[END_REF] due to their design in the empirical risk minimization problem, but evaluating their local calibration conditioned on predicted score values poses challenges [START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF]. In addition, when using more opaque models, such as Random Forest (RF) or Neural Networks, the interpretability of calibration becomes more nuanced, with differing views on their potential (mis)calibration [START_REF] Niculescu-Mizil | Predicting good probabilities with supervised learning[END_REF][START_REF] Guo | On Calibration of Modern Neural Networks[END_REF][START_REF] Hänsch | Stacked Random Forests: More Accurate and Better Calibrated[END_REF][START_REF] Krishnan | Improving model calibration with accuracy versus uncertainty optimization[END_REF][START_REF] Minderer | Revisiting the calibration of modern neural networks[END_REF]. Consequently, various metrics and post-processing recalibration methods, including Platt scaling [START_REF] Platt | Probabilistic outputs for support vector machines and comparisons to regularized likelihood methods[END_REF], isotonic regression [START_REF] Zadrozny | Transforming classifier scores into accurate multiclass probability estimates[END_REF], and Beta calibration [START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF], have been proposed to measure and correct poor calibration in classification models, particularly in scenarios requiring confidence scores.

Considering the comprehensive studies on the concept of calibration, determining the most appropriate metric or recalibration method for a specific dataset and its trained algorithm is not straightforward. Various calibration metrics diverge significantly in their evaluation of models, and a consensus has yet to be established, even in binary tasks. In this regard, after presenting numerical and graphical tools to measure the calibration of a binary predictive model in Section 2, we propose simulating a synthetic dataset for which the true distribution of probabilities is known. By deliberately manipulating these probabilities, we create distorted versions to emulate uncalibrated scores that might be produced by an ML model. With access to the true distribution, we precisely measure the miscalibration of the distorted probabilities using the Mean Squared Error (MSE). Consequently, we can identify calibration metrics that closely resemble the MSE, and introduce a novel metric, the Local Calibration Score (LCS), which uses local regression techniques. Subsequently, in Section 3, we present recalibration approaches to address poor calibration of distorted probabilities. Despite the intentional miscalibration of these scores, the traditional performance metrics have not deteriorated.

The analysis of synthetic data yields insights favoring the novel calibration metric LCS. Subsequently, we compute this measure in a real-world scenario discussed in Section 4, in which a Random Forest algorithm is employed to predict the risk of default. We train both an RF classifier and regressor on the dataset, with the latter demonstrating superior accuracy and calibration, contrasting the results reported by [START_REF] Boström | Calibrating random forests[END_REF]. While selecting the RF algorithm and optimizing its hyperparameters using data from Yeh (2016), we recreate a scenario in which decision makers prioritize finely tuned goodness-of-fit metrics. During this procedure, caution is advised to avoid compromising calibration for the sake of discriminative capacity, particularly with the aim of using model predicted scores.

Our contributions can be summarized as follows:

• In the case of binary regression within our Data Generating Process (DGP), enabling exact calibration calculation, the Expected Calibration Error (ECE) does not emerge as the most robust calibration metric, as well as the Brier score.

• Based on visualization techniques that involve local polynomial regression for calibration curves, we introduce a novel calibration metric named LCS, the relevance of which is validated through the assessment of ground-truth miscalibration on synthetic data.

• When observing the progression of the novel calibration metric -LCS-for different AUC levels during the optimization of both RF regressor and classifier algorithms, we highlight that integrating a calibration metric in the optimization process is significant if one intends to utilize the scores predicted by the classifier.

Calibration

Consider a binary variable D that takes the value 1 if an event occurs and 0 otherwise. In this context, the probability of the event depends on individual characteristics, i.e., p i = s(x i ), where, with sample size n > 0, i = 1, . . . , n represents individuals, and x i the characteristics. The goal is to estimate this probability using a model, such as a Generalized Linear Model (GLM) or an ML model such as an RF. These models estimate a score ŝ(x i ) ∈ [0, 1], allowing the classification of observations based on the estimated probability of the event. By setting a probability threshold τ in [0, 1], one can predict the class of each observation: 1 if the event occurs, and 0 otherwise. However, to interprete the score as a probability, it is crucial that the model is well-calibrated. For a binary variable D, a model is well-calibrated when [START_REF] Schervish | A General Method for Comparing Probability Assessors[END_REF])

P(D = 1 | ŝ(x) = p) = p, ∀p ∈ [0, 1] , (1)
which is equivalent to:

E[D | ŝ(x) = p] = p, ∀p ∈ [0, 1] . (2)
For example, if D represents a credit default, with D = 1 indicating default and D = 0 indicating no default, a model predicting credit default is well-calibrated if the values of model's predicted probability ŝ(x) closely match the actual observed probability of default. This means that if the model predicts a probability of 0.8 for a particular individual, the actual default rate for individuals with a predicted probability of 0.8 should be close to 0.8. It should be mentioned that conditioning by {ŝ(x) = p} leads to the concept of (local) calibration; however, as discussed by [START_REF] Bai | Don't just blame overparametrization for over-confidence: Theoretical analysis of calibration in binary classification[END_REF], {ŝ(x) = p} is a.s. a null mass event. Thus, calibration should be understood in the sense that

E[D | ŝ(x) = p] a.s. → p when n → ∞ ,
meaning that, asymptotically, the model is well-calibrated, or locally wellcalibrated in p, for any p. From the dominated convergence theorem, it also signifies that, "on average," the model is well-calibrated. To account for this, we will consider multiple replications of finite samples in the simulation study in Section 2.2.

Measuring Calibration

To assess the calibration of a predictive model, the literature provides both metrics and visual approaches for evaluating ŝ. Given the continuous nature of the score with a null mass event, various methods have been proposed to detect (mis-)calibration. We will briefly introduce the most popular measures before presenting a new calibration metric (refer to the LCS in Section 2.1.2).

Quantile based measures

Calibration curve In the binary case, based on Equation 2, constructing a calibration curve to visualize the calibration of a model involves estimating the function g(•) that measures miscalibration on its predicted scores ŝ(x):

g : [0, 1] → [0, 1] p → g(p) := E[D | ŝ(x) = p] . ( 3
)
The g function for a well-calibrated model is the identity function g(p) = p. In practice, from this real-valued setting, it is challenging to have a sufficient number of observations in the training dataset with identical scores to effectively measure calibration defined in Equation 1, resulting in a lack of robustness in the estimation process of these probabilities. A common method for estimating calibration is to group observations into B bins, defined by the empirical quantiles of predicted values ŝ(x). The average of observed values, denoted db with b ∈ {1, . . . , B}, in each bin b can then be compared with the central value of the bin. Thus, a calibration curve can be constructed by plotting the centers of each bin on the x-axis and the averages of corresponding observations on the y-axis, also referred to as reliability diagrams [START_REF] Wilks | On the combination of forecast probabilities for consecutive precipitation periods[END_REF]. When the model is well-calibrated, all B points lie on the bisector.

Expected Calibration Error Given a sample size n, the Expected Calibration Error (ECE) [START_REF] Pakdaman Naeini | Obtaining well calibrated probabilities using bayesian binning[END_REF] is determined using two metrics within each bin b ∈ {1, . . . , B} of quantile-binned predicted scores ŝ(x): accuracy acc(b), which measures the average of empirical probabilities or fractions of correctly predicted classes, and confidence conf(b), indicating the model's average confidence within bin b by averaging predicted scores. The ECE is then computed as the average over the bins using:

ECE = B b=1 n b n | acc(b) -conf(b) |
where n b is the number of observations in bin b. Given that each bin b is associated with set I b containing the indices of instances within that bin,

acc(b) = 1 n b i∈I b 1 di =d i and conf(b) = 1 n b i∈I b ŝ(x i ) ,
are, respectively, the accuracy and the confidence of the model in bin b. The predicted class di for observation i is determined based on a classification threshold τ ∈ [0, 1], where di = 1 if ŝ(x i ) ≥ τ and 0 otherwise. Notably, the aforementioned ECE corresponds to the exact definition of calibration for multi-class prediction [START_REF] Guo | On Calibration of Modern Neural Networks[END_REF]. For a recent study investigating the application of this measure to geographic data with ordinal classes, we refer to [START_REF] Machado | Geospatial disparities: A case study on real estate prices in paris[END_REF].

Local Regression based measure

We propose an alternative approach to visualize model calibration, aiming for a smoother representation than that provided by the method based on quantiles.

Smoothed calibration curve

Instead of defining bins, we estimate calibration using a local regression. The benefit of employing local regression techniques over bin-based visualization lies in the fact that given the number of data points, the precision of quantile binning can be suboptimal when determining the appropriate bin count. By contrast, with local regression, one can specify the percentage of nearest neighbors, providing greater flexibility.

Local regression involves applying local polynomial regression techniques, with the definition of "local" adopting various forms, including bandwidth or the count of nearest neighbors. These characteristics are determined using different approaches, enhancing the precision and smoothness of the analysis. The degree of the polynomial determines the construction of the polynomial regression. For this purpose, we will use the locfit package from R [START_REF] Loader | Fitting with LOCFIT[END_REF]. The locfit library selects a specific set of evaluation points from the dataset to conduct the regression fit. By default, the evaluation structure follows a tree-based pattern: the algorithm confines the dataset within a rectangle and divides it into two equalsized rectangles. Polynomial parameters from the fit at selected points are then used to interpolate values at other locations. While this approach may have been disregarded in high dimensions due to poor properties, it is highly efficient in small dimensions, as in this case with only one predictive feature, ŝ(x). Moreover, this approach is backed by a solid theoretical foundation, with well-established statistical properties.

To obtain a calibration curve, we fit a local regression model with a degree of 0.1 We perform a regression of observed events on predicted probabilities. Then, we employ the trained model to make predictions on a linear space with values in [0, 1]. This approach provides a continuous calibration profile, based on uniformly distributed points to evaluate g.

The local regression model relies on neighboring points to make predictions. However, problems arise at the edges when the range of predicted scores ŝ(x) does not cover the entire interval [0, 1]. In such cases, predicted values beyond this range may deviate from the bisector, leading to a misinterpretation of calibration. To prevent this issue, we adjust the linear space used for predictions by the local regression model to align with the full range of observed scores ŝ(x).

Local Calibration Score Similar to how the Integrated Calibration Index is defined using the LOESS regression method [START_REF] Austin | The integrated calibration index (ICI) and related metrics for quantifying the calibration of logistic regression models[END_REF], we introduce our methodology, the Local Calibration Score (LCS), which is based on the calibration curve constructed using locfit, as detailed in Section 2.1.2. The calculation of the LCS relies on the disparities between this curve and the bisector, in the range from 0 to 1, weighted by the density of the predicted scores ŝ(x). To execute this, following the methodology outlined in Section 2.1.2, a local regression of degree 0, denoted as ĝ, is fitted to the predicted scores ŝ(x). This fit is then applied to a vector of linearly spaced values within the interval [0, 1]. Each of these points is denoted by l i , where i ∈ {1, . . . , N }, with N being the target number of points on the visualization curve. The LCS is calculated by averaging the squared differences between each predicted score ĝ(l i ) and its corresponding linearly spaced value l i , weighted by the density of the observed scores at l i , corresponding to w i :

LCS = n i=1 w i ĝ(l i ) -l i 2 . (4)

Impact of a Poor Calibration

In this section, we aim to explore the impact of a poorly-calibrated model. 2 To do this, we create synthetic datasets where we know the true probability of a binary event occurrence in advance. Instead of relying on model-derived scores ŝ(x), we apply transformations directly to the probabilities. These transformed values serve as surrogate scores, representing what could be obtained from an ML model. This approach enables us to assess the effects of poor calibration on different calibration metrics, particularly our novel calibration measure LCS, as well as on various calibration curves. Additionally, it allows us to investigate the impact of poor calibration on standard goodness-of-fit metrics.

Synthetic data Following [START_REF] Gutman | Propensity score models are better when post-calibrated[END_REF], we consider a binary variable assumed to follow a Bernoulli distribution:

D i ∼ B(p i )
, where p i is the probability of observing D i = 1. We define p i using the sigmoid function:

p i = 1 1 + exp(-η i ) , (5)
Here, η i is defined by the equation:

η i = a 1 x 1 + a 2 x 2 + a 3 x 3 -a 4 x 4 + ε i , (6)
where x 1 , x 2 , x 3 , and x 4 are randomly drawn from a uniform distribution U [0, 1], where (a 1 , . . . , a 4 ) are scalars arbitrarily set in our case to (a 1 , a 2 , a 3 , a 4 ) = (0.1, 0.05, 0.2, -0.05), and where ε i ∼ N (0, 0.5 2 ). We generate n = 2, 000 observations using this DGP.

Of particular importance, with this experimental setup and a well-defined framework, we establish in the following proposition that logistic regression is asymptotically well-calibrated. The proof is detailed in Appendix B.1. Proposition 2.1. Consider a dataset {(d i , x i )}, where x are k features (k being fixed), so that D|X = x ∼ B s(x) where

s(x) = exp[β 0 + x ⊤ β] 1 + exp[β 0 + x ⊤ β]
.

Let β 0 and β denote maximum likelihood estimators. Then, for any x, the score is defined as

ŝ(x) = exp[ β0 + x ⊤ β] 1 + exp[ β0 + x ⊤ β] is well-calibrated in the sense that E[D | ŝ(x) = p] a.s. → p as n → ∞.
Proof. Proof in Appendix B.1. If k were increasing with n, [START_REF] Bai | Don't just blame overparametrization for over-confidence: Theoretical analysis of calibration in binary classification[END_REF] showed that logistic regression is over-confident. However, assuming here that k is fixed provides a complementary perspective. This justifies the possibility of achieving a perfectly calibrated model with this DGP, making the proposed synthetic data appealing to study (mis-)calibration.

Next, we introduce two types of transformations to the true probabilities p to simulate uncalibrated modeling: one directly applied to the latent probability and another applied to the linear predictor η. Specifically, we introduce a scaling parameter that modifies the latent probability, altering Equation 5 to:

p u i = 1 1 + exp(-η i ) α .
(7)

A second scaling parameter which modifies the linear predictor changes Equation 6 to:

η u i = γ × (-0.1)x 1 + 0.05x 2 + 0.2x 3 -0.05x 4 + ε i , (8)
The resulting values, given by p u , are considered as the scores ŝ(x) that could be returned by a predictive model. Note that when α = γ = 1, no transformation occurs, representing the benchmark situation of a well-calibrated model. In the simulations, we examine variations in α and γ across the range 1/3, 1, 3, considering each parameter individually while keeping the other fixed at 1. The effects of the transformations on the probabilities are shown in Figure 1.3 Notably, decreasing (increasing) α shifts values closer to 1 (to 0). Decreasing γ concentrates values around 0.5, while increasing γ disperses probabilities around 0.5. For each value of parameter α or γ, we generate 200 independent samples, each consisting of 2, 000 observations. Within each sample, after applying p u , we compute the various metrics previously mentioned. Since we are aware of the true scores, we calculate the associated MSE, n -1 n i=1 (p i -p u i ) 2 , representing the ground-truth calibration and enabling us to understand the characteristics that accurate calibration metrics should exhibit. Its equivalent when the true probability distribution is not known corresponds to the Brier score [START_REF] Brier | Verification of forecasts expressed in terms of probability[END_REF], where p i is replaced by d i . Figure 2 presents the results in the form of boxplots for the calculation of the calibration measures when α varies (top) or when γ varies (bottom). For each metric, a degradation is observed when the scores deviate from the true probabilities, i.e., when α = 1 or γ = 1. It is noteworthy that, with data generated from our DGP, the proposed LCS approach remarkably reflects the true error being close to 0 when α = γ = 1. While the Brier score mirrors the dynamics of the MSE, it registers excessively high values, surpassing 0.23, for the true probability distribution. Conversely, the ECE appears less suitable for assessing calibration within this binary scenario. Calibration can also be assessed using calibration curves. Figure 3 shows for each type of score distortion the average of the calibration curve computed using local regression with locfit, for the 200 simulations. The error band corresponds to a 95% bootstrap confidence interval. The distribution of true probabilities is depicted at the top of each graph. 4 The curves for cases where the model is effectively well-calibrated (α = 1 or γ = 1) are notably aligned with the bisector. Conversely, in other scenarios, the calibration curves deviate considerably from the alignment with the bisector. Although similar visualizations are obtained with the calibration curves computed using bins defined by quantiles (see Figure 10 in the Appendix), calibration curves obtained with local regression exhibit smoother patterns.

Finally, we explore the sensitivity of standard goodness-of-fit metrics to calibration. The applied transformations to probabilities for defining scores reflecting a poorly calibrated model show no discernible impact on standard goodness-of-fit measures, as illustrated in Figure 4. Metrics such as accuracy, sensitivity, specificity, or AUC are not impacted by monotone transformations on the probabilities. When the focus shifts to using predicted scores from a binary classifier rather than merely class prediction accuracy, standard goodness-of-fit metrics may not detect poor calibration. This section explores the benefits of the proposed DGP, highlighting the rationale for employing a post-hoc recalibrator to address this issue.

In the following section, we examine different recalibration approaches, drawing from various techniques in the literature. We then delve deeper into our understanding of local regression methodologies for recalibration in Section 3.1.1.

Recalibration

When using scores generated by a model predicting the probability of a binary event, the literature advocates recalibrating the model by applying a transformation to the scores [START_REF] Platt | Probabilistic outputs for support vector machines and comparisons to regularized likelihood methods[END_REF][START_REF] Zadrozny | Transforming classifier scores into accurate multiclass probability estimates[END_REF][START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF].

To address overfitting during the learning of this mapping, it is recommended to partition data into three sets [START_REF] Zadrozny | Transforming classifier scores into accurate multiclass probability estimates[END_REF][START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF]: a training set for classifier training, a calibration set for recalibrator training, and a test set for computing calibration metrics. Adequate data is required for this tripartite division.

In this section, we present the prevalent recalibration techniques found in the literature. Then, we evaluate their impact on calibration performance using metrics and visualization methods from Section 2.

Recalibration Methods

Platt Scaling This approach was initially introduced to map SVM outputs to well-calibrated posterior probabilities [START_REF] Platt | Probabilistic outputs for support vector machines and comparisons to regularized likelihood methods[END_REF]. The method consists of applying Logistic Regression to the output probabilities of a binary classifier. The Logistic Regression parameters, a and b, are learned on a calibration dataset, held out from a train dataset:

g ŝ(x) = 1 1 + exp -aŝ(x) + b . ( 9
)
Isotonic Regression This solution arises from a constrained optimization problem [START_REF] Zadrozny | Transforming classifier scores into accurate multiclass probability estimates[END_REF], solved using the Pool-Adjacent-Violators Algorithm, ensuring that corrected predicted scores remain monotonic. In Equation 10,

d (i) corresponds to the value in {d 1 , • • • , d n } associated with the i-th largest predicted score {ŝ(x 1 ), • • • , ŝ(x n )}. min β 1 ,...,βn n i=1 (d (i) -β i ) 2 s.t. β 1 ≤ . . . ≤ β n . ( 10
)
Isotonic regression assumes the initial model's predicted scores ŝ(x) are wellordered, limiting its ability to correct non-monotonic probability distortions, as seen in methods like Platt scaling (Equation 9).

Beta Calibration

The Beta calibration method [START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF] estimates the following calibration curve using three parameters, a, b and c:

g ŝ(x) = 1 1 + exp {-c} ŝ(x) a /(1 -ŝ(x)) b . ( 11
)
The condition a, b ≥ 0 leads to an increasing map function. In contrast to Platt scaling, the Beta calibration family includes the identity function, allowing it to maintain score calibration when it is already calibrated.

Local Regression

Interestingly, the local regression method serves a dual role, functioning both as a visualization tool (as discussed in Section 2.1.2) and as an approach for recalibration. Specifically, when fitting a local regression of degree 0 to the predicted scores ŝ(x), it estimates the following quantity:

Ê[D | ŝ(x) ∈ V p ],
where p ∈ [0, 1]. Here, V p represents the neighborhood of a given p, defined using a percentage of nearest neighbors among the set of evaluation points when using locfit (as mentioned in Section 2.1.2). By definition of local regression of degree 0, the estimation of the expected value approximates the function g defined in Equation 2. Nevertheless, to improve the smoothness of expectancy estimations, particularly in regions with limited data points (e.g., when the predicted scores ŝ(x) are close to 0 and 1), employing locfit with polynomial degrees of 1 and 2 is an alternative worth considering, while remaining aware of the potential occurrence of negative values.

Scores Recalibration

We apply the recalibration techniques presented in Section 3.1 to our simulated datasets. For both calibration and test samples, we calculate calibration and goodness-of-fit metrics after replacing the uncalibrated scores p u with their recalibrated values, denoted p c . For each metric and scenario where we transformed the true probabilities to induce poor calibration, we compute the difference between the metric obtained with the recalibrated scores p c and that obtained with the uncalibrated scores p u . The results for MSE, AUC and LCS are presented in Figure 5, specifically for γ = 3. 5Regardless of the method employed, recalibration leads to a decrease in MSE in both the calibration and test samples. However, this appear to occur at the expense of a reduction in AUC in the test sample. Meanwhile, calibration, measured with the LCS, indeed improves after applying any of the considered techniques. The same observation holds true for the calibration curves calculated using the method based on local regression. 6Overall, we observe that across our 200 simulations, the use of recalibration techniques allows for achieving improved calibration, which comes with sacrificing a slight decrease of precision. However, caution must be exercised when interpreting these results. The subsequent section seeks to extend the analysis by using a real-world dataset, thereby verifying whether the observed outcomes in simulations are generalizable and not contingent on the specifics of the DGP.

Experimentation : Calibration of a Random Forest

In this section, we focus on the application of RF models to predict the occurrence probability of a binary event D ∈ {0, 1}. First, we compare RF when considered in the context of regression versus calibration. Then, we investigate the relationship between seeking accurate class predictions and the calibration of models.

Classifier or Regressor

When RFs are employed to predict binary outcome variables, one can either train a classifier or a regressor. However, the scores returned by these two types of models are computed in very different ways. The score returned by an RF classifier corresponds to a majority class vote, whereas an RF regressor estimates probabilities by averaging initial class membership probabilities.

In a classification forest comprising M trees, each tree m ∈ {1, . . . , M } of the forest returns, for an observation i, a majority vote dm,i ∈ {0, 1}. The predicted score for this observation is the average of the majority votes within the tree, pclass i = M m=1 dm,i . When training an RF regression on a binary target variable, the score returned by the forest for an observation i is preg i = M m=1 pm,i , where pm,i is calculated as the average of the observations in the terminal leaf to which observation i belongs to the tree m.

Given the theoretical reliance of the forest regressor on probabilities, we anticipate that its calibration should be superior to that of the classifier, compared to [START_REF] Boström | Calibrating random forests[END_REF]. Our objective is to compare the goodness-of-fit and calibration metrics of both models, both with and without the application of recalibration techniques on the estimated scores. Data For our illustrations, we use data obtained from UCI (Yeh, 2016), presenting research customers' default payments in Taiwan. This dataset contains n = 30, 000 instances and 23 numeric features. The outcome variable, corresponding to the observed default payment in next month, is positive in 22.12% of cases. Following the methodology outlined in [START_REF] Subasi | Prediction of default payment of credit card clients using data mining techniques[END_REF], we employ the Synthetic Minority Over-sampling Technique (SMOTE) [START_REF] Chawla | Smote: synthetic minority over-sampling technique[END_REF] at a rate of 200% to rebalance the data.

To predict default payments, we apply both RF algorithms on the rebalanced dataset. 7 The training set comprises 50% of the dataset, while the remaining portion is further divided into a calibration set used to train a recalibrator and a test set to evaluate the goodness-of-fit and calibration of the models. Given the size of the dataset, this approach is feasible. For smaller datasets, Park and Ho (2020) recommend using out-of-bag samples for training the recalibration method.

Methodology

We conduct a grid search to find the set of hyperparameters that optimize a criterion: the out-of-bag MSE for the regressor, and the error rate for the classifier. The hyperparameters we vary include the number of trees, the number of variables considered for splitting, and the minimum number of observations in terminal nodes. 8 Once the hyperparameters are selected for both types of forests, the forests are used to make predictions on the remaining data not used in training. We then perform 200 simulations by splitting these predictions into two samples: a calibration sample and a test sample. For each simulation, we measure calibration and goodness-of-fit before and after recalibration, on both samples. This process yields a distribution of estimation quality and calibration metrics.

Calibration Simulation results are presented in Figure 6. Considering the inefficiency of the ECE and Brier score, observed in the previous Section (see Figure 2), our novel metric, LCS (Section 2.1) is the only calibration metric shown in the Figure . For the goodness-of-fit metric, we only calculate the AUC, as it is invariant to the probability threshold.9 Also, the true MSE cannot be computed here, due to the unavailability of true data probability distributions.

Figure 6 illustrates that prior to recalibration, the regressor exhibits slightly higher AUC and superior calibration than the classifier, thereby confirming the hypothesis outlined in Section 4.1. Moreover, although both algorithms appear to be well-calibrated, all recalibration methods enhance the calibration of these models without compromising the AUC. However, it is worth noting that isotonic regression seems to overfit in terms of AUC on the test set, an observation also documented by [START_REF] Kull | Beyond sigmoids: How to obtain well-calibrated probabilities from binary classifiers with beta calibration[END_REF].

Optimizing Class Predictions and Calibration

ML models are typically fine-tuned to optimize hyperparameters to maximize goodness-of-fit in binary classification, often with less emphasis on calibration. However, the models that achieve the highest accuracy may not necessarily exhibit superior calibration. To investigate this issue, we assess the LCS across various goodness-of-fit levels, as measured by AUC, throughout the hyperparameter optimization process. Figure 7 illustrates that for both types of RF, the order of calibration and performance metrics on the train sample is respected by the test sample. Furthermore, this Figure demonstrates that, optimizing the performance of an RF classifier with respect to the AUC reduces the calibration on the test sample. Thus, when using this type of model to directly employ its predicted scores, it may be necessary to consider both a performance metric and a calibration measure in the hyperparameter optimization process to consider these scores as the probabilities of belonging (right). Each point represents an estimation obtained from a set of hyperparameters. The gray lines help identify the estimations made within each sample using the same model.

to the positive class. In contrast, in the case of regression, optimizing the AUC appears to enhance calibration on the test sample. However, opting for a model with an AUC not too close to 1 could be more advantageous for preventing overfitting. This involves sacrificing a small percentage of the AUC on the test sample while preserving good calibration, as depicted by the test points falling within an AUC range of 0.89 to 0.92.

Conclusion

This study aims to deepen our understanding of various calibration measures and methods for recalibrating binary classifiers. This is achieved by analyzing a simulated dataset generated from a logistic function, where the true probability distribution is known. We highlight the flexibility of synthetic data, unveiling nuances in calibration metrics, thereby identifying limitations in the Brier score and ECE within this context. We correct these limitations by introducing a novel calibration metric: the Local Calibration Score. This underscores the importance of local regression techniques for visualization and classifier recalibration. Experimental results using RF classifier and regressor to predict default risk demonstrate slightly better accuracy and calibration with the regressor. Notably, our evaluation of the LCS across various AUC levels during hyperparameter optimization reveals that RF classifiers achieving the highest goodness-of-fit do not necessarily exhibit superior calibration.

From Uncertainty to Precision: Enhancing Binary Classifier Performance through Calibration (Supplementary Material)

A Online Replication Book

We provide replication material on GitHub https://github.com/fer-agathe/calibration_binary_classifier.

B Proofs

B.1 Proof of Proposition 2.1

As a starting point, suppose that there is single feature, so that x can be denoted x.

Suppose that D|X = x ∼ B s(x) where

s(x) = exp[β 0 + β 1 x] 1 + exp[β 0 + β 1 x]
with β 1 = 0, and let β 0 and β 1 , denote maximum likelihood estimators, so that the model is well specified. Then, for any x, the score is

ŝ(x) = exp[ β0 + β1 x] 1 + exp[ β0 + β1 x] .
Here both s and ŝ are continuous and invertible. Given p ∈ (0, 1), {ŝ(x) = p} corresponds to {x = ŝ-1 (p)}, since mapping ŝ is one-to-one. Thus, D conditional on {ŝ(x) = p} is therefore a Bernoulli variable with mean s(ŝ -1 (p)). And because ŝ and s are continuous and bijective functions

β0 → β 0 β1 → β 1 as n → ∞ =⇒ ∀x, p ŝ(x) → s(x) ŝ-1 (p) → s -1 (p) as n → ∞
since the model is well specified, where the convergence is in probability here. Thus

s(ŝ -1 (p)) → p as n → ∞ and E[D|ŝ(x) = y] = p(ŝ -1 (p)) → p as n → ∞.
This property holds in higher (fixed) dimensions, if the model is well specified, since functions are continuous and invertible (linear models with continuous and invertible link functions). The case where the dimension of x increases with n is discussed in [START_REF] Bai | Don't just blame overparametrization for over-confidence: Theoretical analysis of calibration in binary classification[END_REF].

C Calibration with Simulated Data

This appendix provides additional figures regarding the 200 simulations made using the data generating process described in the main text.

We generated data from the data generating process described in Equations 5 and 6 and applied transformations either on the probabilities (varying α in Equation 7) or on the linear predictor (varying γ in Equation 8).

For each value of α = {1/3, 1, 3} and γ = {1/3, 1, 3}, we generated 200 datasets. The distribution of a single dataset for each value of α and γ is shown in Figure 8. In each case, the metric is computed using either the true probabilities p which represents the ground truth, the uncalibrated scores p u , or the recalibrated scores, p c . The recalibration techniques are the following: Platt scaling, Isotonic regression, Beta calibration, and Local regression with varying degrees: degree 0, degree 1, and degree 2. See Table 1 for a summary. This figure complements Figure 5 from the main text. However, unlike in Figure 5, the values are not expresses here as the difference with the value observed in the case where the scores used are the uncalibrated scores p u . The calibration curves computed after recalibrating the scores are calculated on the calibration set (shown in orange), and on the test set (shown in green), for each poor calibration scenario, where α or γ vary. When α = 1 or γ = 1, the transformed scores p u are in fact equal to the true probabilities p. We consider different scenarios. In Figure 11a, the scores are the true probabilities p and no recalibration technique is employed. In Figure 11b, the scores are the uncalibrated values p u and no recalibration technique is employed either. Then, we consider recalibration techniques on these uncalibrated values: Platt scaling (Figure 12a), isotonic regression (Figure 12b), beta calibration (Figure 12c), and local regression (Figure 13a for degree 0, Figure 13b for degree 1, and Figure 13c for degree 2).

C.1 Metrics

C.3 Calibration Curves with Local Regression

In addition to the calibration curves computed using the quantile binning approach, we provide the calibration curves computed on the simulations using local regression on the calibration set (shown in orange) and on the test set (shown in green), depending on the scores used to obtain the calibration curves: true probabilities (Figure 14a), uncalibrated transformed probabilities (Figure 14b), transformed probabilities recalibrated with Platt Scaling (Figure 15a), transformed probabilities recalibrated with isotonic regression (Figure 15b), transformed probabilities recalibrated with beta calibration (Figure 15c, and transformed probabilities recalibrated with local regression with degree 0 (Figure 16a), with and with degree 2 (Figure 16c).

D Calibration of a Random Forest

We performed a grid search to select the best set of hyperparameters to train two types of forests on a train set with 50% of the observations: regression forests, and classification forest. For each type of forest, we varied the following hyperparameters:

• ntree, the number of trees: 100, 300, 500

• mtry, the number of variables to consider for a split: 1, 2, . . . , 12, where 12 represents half the number of features

• nodesize, the minimum size in terminal nodes: 5, 10, 15, 20.

In total, this corresponds to training 144 regression forests and 144 classification forests. For the regression forest, we computed the out-of-bag MSE, whereas for the classification forest we computed the out-of-bag error rate (computed as the number of incorrectly classified observations over the total number of observations). The best set of hyperparameters was selected as the one which minimizes the out-of-bag criterion. For the regression forest, the best set of hyperparameters turned out to be ntree=500, mtry=5, and nodesize=5. For the classification forest, the best set was ntree=300, mtry=5, and nodesize=5.

After obtaining these hyperparameters, we split the remaining 50% of observations into a calibration and a test sets of equal size. This split was performed randomly over 200 different replications. In each replication, we trained a recalibrator in the calibration set, using the predicted scores from the forest. As with the simulated data, we considered the following recalibration techniques: Platt scaling, Isotonic regression, Beta calibration, and Local regression with varying degrees: degree 0, degree 1, and degree 2. We then computed goodness-of-fit metrics and calibration metrics on the calibration and on the test set, for each forest and each replication, using either the uncalibrated scores ŝ(x) (denoted as p u ) or on the scores recalibrated with the different methods.

D.1 Metrics

The different metrics are shown in Figure 17 for goodness-of-fit, and in Figure 18 for calibration. Both Figures complement Figure 6 from the main text.

Figure 1 :

 1 Figure 1: Distorted Probabilities as a Function of True Probabilities, Depending on the Value of α (left) or γ (right).

Figure 2 :

 2 Figure 2: Calibration Metrics on 200 Simulations for each Value of α (top) or γ (bottom).

Figure 3 :

 3 Figure 3: Calibration Curve Obtained with Local Regression, on 200 simulations for each Value of α (top) or γ (bottom). Distribution of the true probabilities are shown in the histograms (gold for d = 1, purple for d = 0).

Figure 4 :

 4 Figure 4: Standard Goodness of Fit Metrics on 200 Simulations for each Value of α (top) or γ (bottom). The probability threshold is set to τ = 0.5.

Figure 5 :

 5 Figure 5: Metrics After Recalibration (for γ = 3), on the Calibration (transparent colors) and on the Test Set (full colors).

Figure 8 :

 8 Figure 8: Distribution of Distorted Probabilities p u Depending on α (top) or γ (bottom).

Figure 9

 9 Figure9displays different calibration metrics (True MSE in Figure9a, the Brier Score in Figure9b, the Expected Calibration Error in Figure9c, and the Local Calibration Score in Figure9d), for the complete set of simulations. Each graph of each panel in the Figure shows the distribution of a metric computed over 200 replications of the simulations for a scenario in which the scores are distorted, on the calibration set (transparent colors) and on the test set (solid colors). In each case, the metric is computed using either the true probabilities p which represents the ground truth, the uncalibrated scores p u , or the recalibrated scores, p c . The recalibration techniques are the following: Platt scaling, Isotonic regression, Beta calibration, and Local regression with varying degrees: degree 0, degree 1, and degree 2. See Table1for a summary. This figure complements Figure5from the main text. However, unlike in Figure5, the values are not expresses here as the difference with the value observed in the case where the scores used are the uncalibrated scores p u .

Figure 10

 10 Figure 10 illustrates, mirroring Figure 3 the overlay of 200 calibration curves computed using bins defined by quantiles, for each type of score distortion (varying α or γ). The distribution of true probabilities is depicted at the top of each graph of the figure.The calibration curves computed after recalibrating the scores are calculated on the calibration set (shown in orange), and on the test set (shown in green), for each poor calibration scenario, where α or γ vary. When α = 1 or γ = 1, the transformed scores p u are in fact equal to the true probabilities p. We consider different scenarios. In Figure11a, the scores are the true probabilities p and no recalibration technique is employed. In Figure11b, the scores are the uncalibrated values p u and no recalibration technique is employed either. Then, we consider recalibration techniques on these uncalibrated values: Platt scaling (Figure12a), isotonic regression (Figure12b), beta calibration (Figure12c), and local regression (Figure13afor degree 0, Figure13bfor degree 1, and Figure13cfor degree 2).

  Curves Defined with Bins, on 200 Simulations for each Value of α (top) or γ (bottom).

Figure 11 :

 11 Figure 11: Calibration Curves Calculated with True Probabilities (Panel a) or with Uncalibrated Scores (Panel B) as the Scores. The curves are obtained with quantile binning, for the calibration set (orange) and for the test set (green) for varying values of α and γ. The curves of the 200 replications of the simulations are superimposed. The histogram on top of each graph show the distribution of the true probabilities (Panel a) and of the uncalibrated scores (Panel b) (in gray).

Figure 12 :

 12 Figure 12: Calibration Curves Calculated with Scores Recalibrated Using Platt Scaling (Panel a), Isotonic Regression (Panel b), or Beta Calibration (Panel c). The curves are obtained with quantile binning, for the calibration set (orange) and for the test set (green) for varying values of α and γ. The curves of the 200 replications of the simulations are superimposed. The histogram on top of each graph show the distribution of the uncalibrated scores (gray), and that of the calibrated scores (blue).

Figure 13 :

 13 Figure 13: Calibration Curves Calculated with Scores Recalibrated Using Local Regression with Degree 0 (Panel a), Degree 1 (Panel b), or with Degree 2 (Panel c). The curves are obtained with quantile binning, for the calibration set (orange) and for the test set (green) for varying values of α and γ. The curves of the 200 replications of the simulations are superimposed.The histogram on top of each graph show the distribution of the uncalibrated scores, and that of the calibrated scores.

Table 1 :

 1 Recalibration Techniques Used in the Simulations

	Legend	Scores Used	Recalibration Technique
	True Prob.	Simulated true prob. No recalibration technique
	No Calibration	Transformed prob.	No recalibration technique
	Platt Scaling	Transformed prob.	Platt Scaling
	Isotonic	Transformed prob.	Isotonic regression
	Beta	Transformed prob.	Beta calibration
	Locfit (deg = 0) Transformed prob.	Local reg. with degree 0
	Locfit (deg = 1) Transformed prob.	Local reg. 1
	Locfit (deg = 2) Transformed prob.	Local reg. with degree 2
	C.2 Calibration Curves with Quantile Binning

The choice of degree 0 allows us to compute the local mean of observed events in the vicinity of a predicted probability, illustrating calibration from Equation

2 Replication material: https://github.com/fer-agathe/calibration_binary_classifier.

FernandesMachado, Charpentier, Flachaire, Gallic, Hu

See Figure8in the Appendix for histograms.

This is represented as a histogram obtained by averaging the number of observations in each bin across the 200 simulations.

Further values for α and γ are in the Appendix, Figure9.

Calibration curves for all values of γ and α are provided in the Appendices, in Figures 11, 12, and 13 for curves obtained by quantile binning, and in Figures 14, 15, and 16 for curves obtained by local regression.

We used the the randomForest() function from the R package randomForest.

For further details on the grid search, refer to Section D in the appendices.

Additional metrics are reported in the Appendix, in Figure17for goodness-of-fit, and in Figure18for calibration.

Philipp Ratz for his valuable comments. Emmanuel Flachaire and Ewen Gallic acknowledge that the project leading to this publication has received funding from the French government under the "France 2030" investment plan managed by the French National Research Agency (reference: ANR-17-EURE-0020) and from Excellence Initiative of Aix-Marseille University -A*MIDEX. François Hu acknowledges that the project is funded by the Natural Sciences and Engineering Research Council of Canada (NSERC) Emerging Infectious Diseases Modelling Initiative (EIDM), awarded to the Mathematics for Public Health (MfPH) program.


Florent Blancho

Mélanie Davranche

Adrien Léon

Rémi Marsac

Stéphanie Reynauld

Bruno Grassl

Julien Gigault

Mechanistic description of Lead sorption onto nanoplastics

Keywords: Nanoplastics, Metals, Surface, Adsorption, Modeling Environmental Science: Nano

Although it is currently recognized that nanoplastics (NPs) are present in the environment, their ability to carry metals or other contaminants and the processes involved are still poorly investigated.

Yet, one of the significant threats of plastic debris is found in these associated species. The lack of relevant data is directly explained by our inability to collect enough nanoplastics from the environment for experimental investigation. The environmental NPs models that we have recently proposed offer new opportunities to characterize their interactions with metals. We have developed an adsorption model based on lead (Pb 2+ ) adsorption experiments (adsorption isotherm at pH 5.5 for Pb(II)/NPs ratios between 0.5 and 100 mg gC -1 and a pH adsorption edge between pH 2 and 8 for Pb(II)/NPs at 2.5 and 25 mg gC -1 relative to the studied NPs, Ionic strength = 5 mmol L -1 of NaNO 3 ) using various NPs models to describe surface complexation processes and to predict the sorption capacity of NPs. Model hypotheses have suggested that Pb is adsorbed onto NPs via carboxylic sites mainly as mono-ligand and bi-ligand complexes. The empirical model and linear free energy relationship validated the extrapolated stability constants and binding hypothesis. The adsorption capacities of nanoplastics are highly dependent on their surface oxidation state, shape heterogeneity

Environmenta Significance Statement:

The findings presented in this study shed light on several crucial aspects of nanoplastics that underscore their environmental significance. Nanoplastics exhibit a remarkably high surface reactivity with respect to adsorption. This heightened reactivity implies that they can readily interact with various environmental constituents, potentially affecting ecosystems and organisms. Furthermore, the study highlights the distinct metal adsorption capacities between nanoplastics and polymer nanobeads, underlining the importance of understanding the unique characteristics of different nanoplastic materials, such as size and shape, in the context of metal pollution. Additionally, the formation of mono and bi-ligand complexes with lead (Pb) and nanoplastics bearing -COOH groups elucidates the intricate chemical interactions occurring at the nanoplastic-water interface, emphasizing their potential role as carriers or modifiers of metal pollutants in aquatic systems.

Introduction

Plastic pollution is one of the significant challenges of the Anthropocene Epoch. Plastics benefit society in innumerable ways, and their production has increased considerably over the last century.

This results in the significant generation of plastic debris all throughout the world. Within the size continuum of plastic debris, it is now admitted that nanoscale plastics, i.e. nanoplastics (NPs), are present in our various environments [START_REF] Ter Halle | [END_REF][2][3] . Nanoplastics are colloidal particles < 1000 nm presenting a Brownian motion in aqueous systems, with open and heterogeneous shapes and non-homogeneous charged surfaces 4,5 . Due to these colloidal properties, nanoplastics are therefore potentially highly reactive, especially with regards to their surrounding dissolved contaminants, and they raise several urgent questions about their potential impact on these associated species.

As demonstrated for microplastics (MPs) by Brennecke et al. (2016a), Cu and Zn released from antifouling paint were adsorbed by polystyrene (PS) pellets and in higher concentrations by polyvinyl chloride (PVC) pellets 7 . Turner and Holmes (2015) demonstrated that MPs' metals adsorption is driven by the pH and that altered polymer pellets collected on beaches have a higher metal affinity than virgin pellets 8 . They explained these results as being due to the alteration of the surface that produced the functional sites as well as the presence of organic materials. Tang et al. (2021) used nylon MPs collected in the environment to show that Ni, Cu and Zn adsorption onto MPs is controlled by surface complexation and intraparticle diffusion processes. They also highlighted the importance of the O-containing group in metal adsorption. Recently, artificial neural network calculations have suggested that aged MPs have a higher affinity for metals than virgin MPs and that the sorption 3 capacity depends on the initial metal concentrations 10 . Nonetheless, considering the estimation by Besseling et al. (2019) that nanoparticles (i.e., those with sizes below 100nm in this context) could have concentrations 10 14 times greater than microplastics, it follows that nanoparticles, due to their heightened surface area specificity, are expected to play a more significant influence on the speciation, transportation, and impact of metals within diverse environmental ecosystems. 4,12 .

It is well accepted that the interaction between metals and plastic debris, from microplastics to nanoplastics, is a crucial process in their (eco)toxicological risk assessment. When NPs interact with a metal, this could also increase their bioavailability. During the bio-assimilation of NPs, the associated micropollutants could be released and accumulated in living organisms under favorable physicochemical conditions 13,14 . However, the mechanisms that control the interaction of NPs with metals are still unknown. This is partly explained by the need (and difficulty) to collect a significant amount of NPs from the environment for experimental investigation. To meet this challenge, several authors have produced NPs models. Three production methodologies are proposed: i.

Controlled emulsion-polymerization of spherical nanoparticles with or without a functionalized surface 15,16 ii. Mechanical degradation of pristine polymer pellets, resulting in NPs with heterogeneous shapes and functional surface sites 17,18 iii. Mechanical abrasion of highly altered plastics debris sampled in the environment, resulting in NPs with heterogeneous shapes and functional surface sites 19 .

Although spherical nanoparticles are useful for instrumental calibration due to their controlled size, shape, composition and binding site density 15,16 , these models are not environmentally relevant in terms of shape and size polydispersity 4 . NPs produced by mechanical degradation (polymer pellets and altered plastics debris) are polydisperse in size with heterogeneous and open shapes, as observed for environmental NPs 5,18,19 . However, there are essential differences between NPs produced from polymer pellets and plastic debris since they are not submitted to extreme surface photooxidation in the environment, which strongly influences their site density and, therefore, their reactivity. A first 21 showed that the presence of Pb localized on the altered surface layers of microplastics is attributed to its use as an additive during plastic manufacturing rather than originating from the environment. Lead (Pb) was one of the most commonly used metallic additives in plastic formulations.

This work aims to characterize and quantify the adsorption mechanisms between lead (Pb) and various NPs models: carboxylated polystyrene Latex nanoparticles (PSL surfactant ), soap-free carboxylated polystyrene Latex (PSL free ) 15 , polymorph and polysized nanoparticles (m-NPs) 22 , and plastic nanoparticles produced from the altered layer of plastics debris collected in the environment (e-NPs) 19 . To achieve this, we generated experimental datasets through adsorption experiments conducted under various physicochemical conditions such as concentration and pH. These datasets were employed to assess different binding scenarios and to extrapolate sorption parameters. The model's validity was subsequently confirmed by comparing its predictions with experimental datasets found in the literature. In this context, geochemical modeling serves as a predictive tool for estimating speciation and testing mechanistic hypotheses.

Materials and methods

Chemicals and nanoplastics

All used chemicals were of analytical grade. The solutions were prepared with ultrapure water 108 (Milli-Q system, Millipore). The containers used were (i) trace-metal cleaned with 10% (v/v) HNO3

for 24 h at 45°C, (ii) soaked in and then rinsed with ultrapure water for 24 h at 45°C, and (iii) finally dried at 30°C. All solutions were prepared with ultrapure water (18. iii. polymorph and polysized nanoparticles produced from the milling of polystyrene pellets 18 , (i.e. m-NPs); iv. plastic nanoparticles produced from the altered layer of plastic debris collected in the environment 19 (i.e. e-NPs).

Nanoparticle Properties

The hydrodynamic diameter (D h ) of the NPs was determined using Dynamic Light Scattering (DLS) (VASCO Flex, Cordouan Technologies) for suspensions at 10 ppm without Ionic Strength (IS) (Fig. 1). The shape of the NPs was examined using transmission electronic microscopy (TEM) (JEM 2100 HR, Jeol). Transmission electronic microscopy was operated at an acceleration voltage of 200-kV with LaB 6 as an electron source, the point and line resolution were 2.3 Å and 1.4 Å, respectively.

Particles were photographed with a Gatan Orius SC200D camera, and elemental analysis was performed using an EDX Oxford X-Max 80T detector. The TEM observations showed anisotropic e- and m-NPs (Fig. 1C), in agreement with the previous studies 19,22 . The specific surface area (SSA) of the NPs was determined using the BET method using a Gemini VII instrument (Micromeritics) and according the ISO 9277-2022 protocol 23 . For the polystyrene latex particles, the specific surface area was extrapolated from the particle diameter and assuming that the surface is perfectly smooth. Their diameters were estimated through TEM images treated with ImageJ (n=50).

The surface charge was assessed by potentiometric titration using a Titrano unit controlled by the Tiamo software (Metrohm). All solutions were titrated with 0.01 M NaOH (Honeywell Fluka) under a N 2 flux. Before titration, the pH probe was calibrated on a concentration scale using 1mM HNO 3 .

Titrations were performed at 10, 100, and 1000 mM of NaNO 3 and were started below a pH of 4 by adding a fixed volume of 100 mM HNO 3 . The determination of the protonable groups was performed as described by Spadini at al. (2018) (SI 2). The surface charge of the NPs was determined by performing a calculation using PHREEQC implemented by the SIT thermodynamic database (SI 2).

Potentiometric titrations were used to determine the amount of ionizable function relative to the pH on the PSL surface (SI 2).

Adsorption experiments

Adsorption isotherms were carried out for the Pb(II)/NPs ratios between 0.5 and 100 mg gC -1 . The concentration of the NPs was fixed at 50 and 25 mg L -1 for the PSLs, e-NPs, and m-NPs, respectively.

Note that the concentration of the NPs is in gC L -1 . The mass concentration was preferred over the particle number as it can be directly measured. The ionic strength (IS) was fixed with 5 mmol L -1 of NaNO 3 . For adsorption isotherms, once Pb(II) was introduced into the NP suspensions, the pH was initially set and then continuously maintained at a pH of 5.5 ± 0.05. The filter cut-off used was 0.2 µm (filter syringe, PES, Sartorius) for the PSLs and 30 kDa (ultrafiltration device, PES, Sartorius)

for the e-NPs and m-NPs. After filtration, all samples were acidified with 2% HNO 3 for the ICPMS analysis. The pH-adsorption edge experiments were performed with Pb(II)/NPs ratios of 2.5 and 25 mg gC -1 for the PSLs, e-NPs, and m-NPs, respectively. The pH range was investigated from 2 to 7 using HNO 3 or NaOH and the same method described for the isotherm adsorption. The targeted pH The calculated SI, presented in the Fig. S3, are < 0 for all Pb-based solid phases reported. These results suggest no Pb-based solid precipitation for all condition studied. All adsorption experiments were performed in replicates (n=3).

Chemical Analysis

The lead concentrations were measured with a Quadrupole ICP-MS (Agilent Technologies 7700X).

Before quantification of the metals, calibration curves were performed and validated using certified material references (SLRS-6, National Research Council). A rhodium solution was used as an internal standard to correct the instrumental drift and potential matrix effects. The limit of Pb quantification was determined at 0.0038 µg L -1 . Moreover, the samples were diluted if necessary as the upper limit of the Pb(II) calibration curve was 5 µg L -1 . The NPs concentrations were monitored with the carbon concentration using a Total Organic Carbon (TOC) analyzer (TOC-V Shimadzu).

Modeling approach

The study carried out by Blancho et al. 20 used two complementary modeling approaches. These authors demonstrated that REE binding onto PSLs and e-NPs is provided by REE adsorption onto --COOH surface functions developed at the surface of the NPs. They also demonstrated that REE binding occurs through the formation of mono-ligand and multi-ligand complexes related to the surface properties of the NPs. Therefore, we tested these hypotheses using empirical models to evidence the formation of one or two complexes. The hypothesis of the formation of a single surface complex was tested by using a classical Langmuir isotherm model (Eq. 1), while the formation of weak (monoligand complex, Eq. 2) and strong (multi-ligand complex, Eq. 2) surface complexes was tested by using a double-Langmuir isotherm model. The potentially involved mechanisms were tested and parameterized by accounting for electrostatic effects at the water interface of the NPs using a diffuse layer model (DLM). For both models, the modeling was based on the extrapolations of the experimental datasets, and a least-square fit determined all parameters.

The classical Langmuir model considers the possibility of forming one complex with the discrete surface sites of the adsorbent following the equation:

(Eq. 1) [ ] = [ ] 1 + [ ]
Where Q max is the maximum amount of adsorbed Pb (µg g -1 ), and K ads is the Langmuir adsorption constant, [Pb] ads is the equilibrium adsorbed-phase concentration of Pb (µg g -1 ), and [Pb] aq is the equilibrium aqueous-phase concentration of Pb in solution (µg g -1 ).

The double Langmuir model considers the possibility of forming weak and robust complexes with the discrete surface sites of the adsorbent following the equation:

+ (Eq. 2) [ ] = [ ] 1 + [ ] [ ] 1 + [ ]
where Q weak/strong max, Pb is the maximum amount of adsorbed Pb (µg g -1 ) occurring as weak and robust complexes, and K weak/strong is the Langmuir adsorption constant for weak and strong complexes, respectively.

The mechanistic modeling (DLM) was performed using a combination of geochemical (PHREEQC The Fig. 2 presents the proton interaction with the surface of the studied NPs, highlighting of proton reactive sites. As expected, the dissociation of such groups increased with increasing pH. The IS increasing also favored the surface group deprotonation by inducing an increase in the surface potential (X 0 ) 28 . Titrations of both PSLs at 1M IS showed a typical sigmoidal curve [START_REF] Connors | Chemical kinetics: the study of reaction rates in solution[END_REF][START_REF] Chen | [END_REF] . These titration curves indicated that almost all the ionizable groups were dissociated between pH 4 and 8.

The total surface ionizable group, Hs tot , was considered equal to the concentration of the deprotonated group at pH 8. For the PSL surfactant and PSL free , the Hs tot was thus determined to 0.316 ± 2 10 -3 and 0.086 ± 2 10 -3 mmol g -1 , respectively. According to the SSA, the site density is equal to 6.60 and 3.41 sites nm -2 for PSL surfactant and PSL free , respectively. The extrapolation of Hs datasets allowed us to calculate the intrinsic pKa of PSLs (Table 1, for further information please report to the SI 1). The PSL surfactant pKa was estimated to 4.64 which is a coherent value compare to known COOH pKa [START_REF] Martell | Other Organic Ligands[END_REF] . The PSL free pKa was estimated to 4.11, which consistent with acid's pKa, used during its polymerization, estimated to 4.05 [START_REF] Lützenkirchen | [END_REF] . These results are also consistent with the pKa values determined for the carboxylic functional groups of natural organic matter 33 . sites (Fig. 3). However, no plateau was reached, suggesting that site saturation was not reached or that adsorption was controlled by a multilayer or heterogeneous process. Significant adsorption capacity variations were obtained between the PSLs. For [Pb] free = 2ppm, [Pb] ads was 13 and 15 mg g -1 C for PSL surfactant and PSL free , respectively. This difference was consistent with the -COOH site density (Table 1). For the PSLs, the Pb(II) adsorption increased with the pH. Less than 5% of [Pb] tot was adsorbed at pH [ 4 versus 85% at pH Z 7 (Fig. 2). For the e-NPs and m-NPs, at [Pb] free = 0.4 ppm, [Pb] ads was W 7.5 and 2.25 mg g -1 C, respectively (Fig. 3). These adsorbed concentrations were consistent with the total ionizable site density, Hs tot (Table 1). The adsorption was also strongly dependent on the pH. Less than 15% of [Pb] tot was adsorbed for pH [ 4 versus 90% at pH 7. The pHdependence of the adsorbed concentration highlights the competition between H + and Pb(II) for their adsorption onto the surfaces of the nanoparticles (Fig. 4) (red dotted curve). The formation of both complexes was in line with the e-NPs site density, which was the highest among the studied NPs (Table 1). With regards to the results obtained from the empirical models and the binding hypothesis suggested by Blancho et al. (2022), a one-complex model (homogeneous adsorption process) was tested for both PSLs and the m-NPs; a two-complex (weak mono-ligand and strong bi-ligand complexes, heterogeneous adsorption process) model was tested for the e-NPs. Following Blancho et al. 20 , the -COOH functional groups were considered to be the only ones able to bind Pb onto the NPs.

1
Two hypotheses were tested for the one-complex model: the formation of mono-ligand complexes and the formation of bi-ligand models for both PSLs and the m-NPs. The best-fitting model (lowest RMSE) was reported in Fig. 5, except for the m-NPs, for which the model providing the more consistent log K value was retained (SI 5, Fig. S5). A comparison of both mono-and bi-ligand models for both PSLs and the m-NPs can be found in SI 5 and Fig. S5. Table 2 summarizes the fixed and fitted (calculated) binding parameters.

For e-NPs, only the mono + bi-ligand model was tested with regards to the double-Langmuir model result and the previous results of Blancho et al. 20 , who suggested the formation of mono and multi- ligand complexes between the e-NPs and REE. However, with this model, several assumptions had to be made to allow the model to converge (too many parameters to adjust). Therefore, the bi-ligand complex concentration was fixed equal to that determined by the double-Langmuir model, the log K e- NPsCOOPb+ mono-ligand complex was fixed equal to the log K CH3COOPb+ complex (2.38, Smith and Martell, 1989) and the bi-ligand complex log K e-NPs(COO)2Pb was fixed equal to 2*log K mono.

Mechanistic calculations demonstrated that Pb(II) was bound as a bi-ligand complex and a monoligand complex onto PSL surfactant and PSL free , respectively (Fig. 3 and Table 2). These binding hypotheses are in agreement with their site density (Table 1), the Langmuir model calculations, and the previous results 20 obtained for REE. For PSL surfactant , the calculated log K (NPs -COO)2Pb was 3.76, which was close to the log K (CH3 -COO)2Pb = 4.08 [START_REF] Martell | Other Organic Ligands[END_REF] . For PSL free , the calculated log K NPsCOOPb+ was 2.3, in line with the log K CH3COOPb + = 2.68 [START_REF] Martell | Other Organic Ligands[END_REF] . For the m-NPs, Pb(II) was bound as a bi-ligand complex, with log K (CH3COO)2Pb = 4.94. This value is higher than that of (CH 3 COO) 2 Pb with log K = 4.08, suggesting a high affinity for Pb. However, this affinity does not seem to agree with the m-NPs site density, which was the lowest among the studied m-NP S (Table 1).

Nevertheless, although the site density is an important controlling parameter, the geometric distribution of the functional groups on the NPs surface can also influence the formed complex. Among all the studied NPs, m-NPs have the most heterogeneous shape (Fig. 1). Their size was also larger than that of the e-NPs suggesting the highest aggregation of the m-NPs (Table 1 and Fig. 1).

Furthermore, the presence of -COOH functional groups can also be identified within the plastic crackles, potentially leading to a higher local site density. All these factors influence the geometrical distribution of the surface functional group and can explain that the -COOH sites are close enough to favor the formation of the bi-ligand complex with Pb. e-NPs formed mono and bi-ligand complexes with Pb(II), which was in agreement with the Langmuir calculations and binding hypothesis suggested by Blancho et al. 20 for REE, i.e. adsorption onto -COOH surface functions and formation of mono-ligand and multi-ligand complexes. The calculated log K NPsCOOPb+ was 2.72, a value close to the log K CH3COOPb + = 2.68, as referenced in Martell and Smith (1989). Unfortunately, the constraint the e-NPs of Davranche et al. 12 , this parameter was fitted. Modelling calculations determined a -COOH total concentration of 0.114 mmol g -1 , in the same order of magnitude as the e-NPs produced here. The RMSE = 0.23 was slightly higher than for the other NPs but within the same order of magnitude. The model overestimated the bound Pb for the low Pb concentrations (from [Pb] free < 10 -6 mol L -1 ) and overestimated them for acidic pH (pH < 7) (Fig. 5 C andD). As previously explained, NOM was still present in the adsorption experiment carried out by Davranche et al. 12 . This NOM could have influenced the amount of Pb adsorbed by the NPs by complexing Pb and, competing with the NPs, forming aggregates (with and without NPs) at acidic pH. At acidic pH, dissolved OM could have aggregated to form bigger particles retained in the >3KDa fraction and thus bias the adsorbed Pb concentration. The formation of hetero-aggregates could also modify the binding site density of the NPs. However, at high Pb concentrations and pH (Fig. 5 C andD), it was possible to accurately predict [Pb] ads , indicating no hetero-aggregation. Pradel et al. 37 demonstrated that the source of NOM was a significant parameter of the hetero-aggregation. They showed that hetero-aggregates could be formed with alginate but not with charged OM such as humic acid. These calculations show that the model needs to be improved, in particular to take into account the homo and hetero-aggregation that can occur in the presence of organic matter and at acidic pH. But the good reproducibility of the experimental datasets at circumneutral pH validated the binding hypothesis of our modelling approach.

Environmental implications

The model calculations suggested that mono and bi-ligand complexes are formed between Pb and NPs, except for PSL free . For the e-NPs, both mono-and bi-ligand complexes are produced. Both parameters are controlled by the shape heterogeneity and the NPs aggregation (e.g. m-NPs). However, it seems as though the formed complex is controlled by the site density and the geometrical distribution of the -COO -sites. The present results highlight the role of the degree of NPs degradation which depends on the medium and the prevailing physicochemical conditions. In surface water (river and ocean), photooxidation is expected to be the primary process responsible for the alteration of the e-NPs. Photooxidation mainly produces COO -binding sites, but their density will depend on the

Conclusions

To the best of our knowledge, this work is the first to design a surface complexation model to describe and predict the adsorption of Pb by NPs. Following the previous study 12,20 , two binding hypotheses were tested to explain the adsorption of Pb by NPs. For all NPs models, carboxyl complexes dominate lead adsorption, and bi-ligand adsorption occurred for the PSL surfactant , m-NPs, and e-NPs models. The extrapolated log K was consistent with the log K RCOO-Pb and log K (RCOO) 2 Pb + values from the literature, validating the NPs adsorption models (Fig.s 7 and 8). Metal-NPs interactions depend on the binding site density, shape and aggregation. Milled NPs are less reactive to metal than strongly NPs photo-oxidized under environmental conditions. These results also highlighted the influence of the NPs degradation pathway on the NPs shape and binding site density distribution and thus availability to form mono or bi-ligand complexes. Although NPs competition with natural adsorbents is expected to be low and dependent on their concentration, multi-component sorption should be better studied within environmental conditions to define the metals associated with NPs in the environment. Moreover, competition between metals must be studied in further detail to identify metals having the highest affinity for NPs. Such studies would include a mechanistic approach in order to complete the sorption models proposed in this study and to test their reliability in predicting metal sorption in a more complex matrix.

version 2 )Fig. 1 .

 21 Fig. 1. NP size characterization. A) Log-transformed ACF of the e-NPs and m-NPs superimposed on CMRs standards (200, 500 and 900 nm). B) Population size distribution by intensity estimated with an ACF cumulant fit. C-D) TEM pictures of the e-NPs and m-NPs, respectively.

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. View Article Online DOI: 10.1039/D3EN00677H

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. View Article Online DOI: 10.1039/D3EN00677H

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. View Article Online DOI: 10.1039/D3EN00677H 20 release of the metal, it could also favor their transfer during the colloidal phase in ecosystems and living organisms.

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. View Article Online DOI: 10.1039/D3EN00677H

Environmental Science: Nano Accepted Manuscript

	View Article Online
	DOI: 10.1039/D3EN00677H
	1
	2
	3
	4
	5
	6
	7 8 9 Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM.

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. Prior the Pb(II) adsorption on the nanoplastics, the saturation indices (SI) of Pb-based solid phases were calculated for experimental conditions of isotherm and the pH adsorption edge. The calculations were performed using PhreeQC and MINTEQA2 database, excluding any adsorption on NPs but including an equilibrium with the CO 2 from the atmosphere.

	were 2, 3, 4, 4.5, 5, 5.5, 6 and 7.	View Article Online DOI: 10.1039/D3EN00677H
	1
	2
	3
	4
	5
	6 7 8 9	Environmental Science: Nano Accepted Manuscript

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM.

	View Article Online
	DOI: 10.1039/D3EN00677H
	1
	2
	3
	4
	5
	6 7 8 9	Environmental Science: Nano Accepted Manuscript

Table 1 .

 1 Nanoplastic characteristics. Nanoplastic size and surface properties.

	NPs	D h (nm)	PDI	SSA (m² g -1 )	Hs tot (mmol g -1 )	Site density (site nm -2 )	pKa **	RMSE ***
	PSL surfactant	240	0.06	30.4*	0.316 ± 2e -3	6.60	4.64 ± 0.01	3.2e -5
	PSL free	396	0.06	16.0*	0.086 ± 2e -3	3.41	4.11 ± 0.02	1.7e -5
	e-NPs	262	0.22	13.9	Z 0.282 ± 4e -4	Z 13.6	4.76*	-
	m-NPs	352	0.13	28.9	Z 0.037 ± 1e -4	Z 1.19	4.76*	-

* SSA calculated from the particle size determined via TEM pictures (see Fig.

S1

for the diameter size). ** pKa fixed to the acetate pKa

(Martell and Smith 1977)

. *** RMSE based on the pKa extrapolation.

  Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM.

	View Article Online View Article Online View Article Online
	DOI: 10.1039/D3EN00677H
	1 1
	2 2 2
	3 3 3
	4 4 4
	5 5 5
	6 7 8 9 6 7 8 9 6 7 8 9	Environmental Science: Nano Accepted Manuscript Environmental Science: Nano Accepted Manuscript Environmental Science: Nano Accepted Manuscript

DOI: 10.1039/D3EN00677H

Table 2 . Complexation parameters obtained with empirical and mechanistic approaches.

 2

	Langmuir parameters	PSL surfactant	PSL free	m-NPs	e-NPs
	Q (µmol L -1 / mmol g -1 )	4.15 / 0.16	1.15 / 0.04	0.3 / 0.01
	log K	-0.5	0.3	1.8
	RMSE langmuir	0.08	0.02	0.0005
	Q weak (µmol L -1 / mmol g -1 )				7.06 / 0.27
	log K weak				-1.6
	Q strong (µmol L -1 / mmol g -1 )				0.408 /0.016
	log K strong				0.8
	RMSE double-Langmuir				0.0009
	Surface complexation models	PSL surfactant Bi-ligand	PSL free Mono-ligand	m-NPs Bi-ligand	e-NPs Mono-and Bi-ligand
	H stot (mmol g -1 )	0.32 a	0.09 a	0.037 a	0.282 a
	Hs bi-ligand (mmol g -1 )	-	-		0.02 d
	SSA (m 2 g -1 )	30.4 a	16.0 a	28.9 a	13.9 a
	pKa	4.64 b	4.11 b	4.76 b	4.76 b
	Log K mono (Weak)	-	2.3 ± 0.21 c	-	2.72 c
	Log K bi (Strong)	3.76 c	-	4.94 ± 0.06 c	=2*log K mono
	RMSE	0.203	0.03	0.019	0.107

a: fixed from analytical data (Table

1

) b: fixed and determined from potentiometric titrations (Table

1

) c: fitted from the mechanistic model d: fixed to the strong complex density determined with the double-Langmuir model

Environmental Science: Nano Accepted Manuscript Published

  on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM.

	Environmental Science: Nano	Page 16 of 23
		View Article Online
		DOI: 10.1039/D3EN00677H
	1
	2
	3
	4
	5
	6
	7
	8
	9

Published on 12 February 2024. Downloaded by Bibliotheque de L'Universite de Rennes I on 2/19/2024 9:15:15 AM. View Article Online DOI: 10.1039/D3EN00677H

Acknowledgments

We also acknowledge Ludivine Rault for the TEM observations performed at THEMIS 476 (ScanMAT, UMS 2001 CNRS -University of Rennes 1). Dr. Sara Mullin is also acknowledged for proofreading the English content (https://fr.linkedin.com/in/sara-mullin-0b53b112/fr). We would also like to thank Mrs. Fatimah Sulu Gambari and Laurent Lebreton from Ocean Cleanup for providing the altered plastic debris necessary to produce the e-NPs.

residence time of the plastic debris in the ocean and river. In soil, the degradation is very low and occurs mainly through biotic processes [START_REF] Otake | Biodegradation of low)density polyethylene, polystyrene, polyvinyl chloride, and urea formaldehyde resin buried under soil for over 32 years -Otake[END_REF] , with UV penetrating only the very first cm of the uppermost layer of the soil. In this context, it is therefore expected that only a small amount of metal will be bound by e-NPs in the soil. By contrast, in compost made of household wastes and used for soil amendments, the composting process involves an intense alteration/oxidation of the plastics via the production of OH radicals able to produce COO -sites [START_REF] Xing | [END_REF] . In these conditions, NPs can act as metal adsorbents. Nonetheless, nanoparticles (NPs) can serve as effective adsorbents only when their concentration is sufficiently high to compete with other soil adsorbents, such as natural organic matter and Fe-oxyhydroxides. This is particularly applicable in soil lacking an organic horizon or in soils with low iron content. Furthermore, it is essential to take into account the influence of eco-corona formation 40 , which can either modify the availability of NPs' adsorption sites or enhance metal adsorption by creating new binding sites (additive effect).

Aggregation also seems to impact the sorption process at the NPs surface. Intuitively, aggregation is considered a limitation to metal binding in response to the decreasing availability of sites. However, the results suggest that aggregation associated with the shape heterogeneity can favor metal adsorption by controlling the geometrical distribution of the functional group and creating a local high binding site density that favors the formation of multidentate complexes. Nevertheless, one consequence of this aggregation is the enlargement of particle size, leading to a subsequent reduction in the colloidal stability of the nanoparticles (NPs). As a result, while aggregation may facilitate the formation of a stable complex between the metal and NPs, it can also lead to the entrapment of NPs along with the bound metal within soil pores or sediment due to NPs settling.

Lead may also be present in the plastic matrices as additives. It is expected that lead will be released during the photooxidation of the polymer and then partly adsorbed on the altered layer of the MPs that will further be altered as NPs. Catrouillet et al. 21 calculated a median concentration of Pb in MPs close to 0.3 mg Kg -1 C, i.e. 6500 times lower than the Pb potentially adsorbed by the e-NPs. These authors suggested that the Pb adsorbed on the altered surface of MPs originated from Pb used as additives rather than from the environment. Although this adsorption by NPs partly prevents the

Conflicts of interest

There are no conflicts of interest to declare.

Source of financial support

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. This work was supported by