Cédric Garnier
email: cgarnier@univ-tln.fr

Ivanka Pižeta

Stéphane Mounier

Jean Yves Benaïm

Marko Branica

Influence of the type of titration and of data treatment methods on metal complexing parameters determination of single and multi-ligand systems measured by stripping voltammetry

Keywords: Voltammetry, Metal complexing parameters, Data treatment, Titration, PROSECE

A set of simulated experiments was analysed in order to compare the influence of the titration type and of data treatment methods on the accuracy of metal complexing parameters determination for one-and two-ligand systems. The simulated data corresponded to those obtained by anodic stripping voltammetry and were chosen to represent experiments in linear, logarithmic and decade titration modes. The values of preset complexing parameters for one-and two-ligand systems were chosen to fit into the expected experimental range. Random noise was added to the data prior to the treatment. Five different data treatments were applied: Chau-Buffle, Ružić-van den Berg and Scatchard linearisations, and non-linear fitting and PROSECE optimisations.

The investigation has shown that even in the case of a one-ligand system, logarithmic and decade titrations are much better compared to the linear ones. Linearisation methods are in many cases inferior to those using optimisation algorithms. Random noise has a significant influence on the results of linearisation methods as well. For linearisation methods, in the case of a one-ligand system, high correlation has been found for the confidence interval of the calculated parameters and the difference between the preset and the calculated values. This correlation is proposed to be used as an estimation for the results quality in real experiments. PROSECE is by far superior to other methods in most of the cases due to its flexible and powerful mathematical background. It is highly recommended as a tool for data treatment. Construction of "contour-graphs" enables error prediction of the calculated complexing parameters. PROSECE is proposed as an orientation and valorisation tool in real samples analyses.

Introduction

In natural waters, speciation of trace metals is highly linked to the presence of complexing groups, such as natural organic matter [1][START_REF] Buffle | Complexation of Trace Metals in Natural Waters[END_REF][START_REF] Buffle | [END_REF][4]. The determination of their complexing properties, according to a discrete model that provides insight into the complex interactions by determining, e.g. complexing capacities and stability constants (L iT , K i ), has been the subject of many studies for over 30 years [5][6][7].

For this purpose, different methods of data treatment, derived from voltammetric or potentiometric experiments (anodic and cathodic stripping voltammetry: ASV-CSV, ion selective electrode: ISE, . . . ), have been developed. The Earlier the techniques were based on mathematical data transformation [5,[8][9][10][11][12][13][14] aiming to obtain linearisation. Therefore, they are rather easy to apply, especially for a one-metal-one-ligand system. Increase in the precision of analytical techniques and their automation led to two-ligand models. In that case, mathematical transformation of experimental data (Ružić and Scatchard transformations with two ligands [8,10]), followed by partial linearisation of these data, allows the calculation of complexing parameters of ligands. These linearisation techniques are used in the studies aiming to define complexing properties of natural samples regarding trace metals such as Cd, Cu, Pb, etc. [15][16][17][18][19][20][21][22]. However, taking into account the development of high-performance computing tools to define more complex models which integrate various chemical interactions [23][24][25][26][27][28][29][30], data treatment based on optimisation calculation of the complexing parameters by fitting to experimental data has been developed. Using all data points, arbitrariness of dividing a curve into quasi-linear subsections was eliminated [31]. These data treatments use either bulk data or modified data according to the Ružić or Scatchard approaches [32,33]. Furthermore, they allow the determination of more than two different groups of ligands [34], which seems to be more useful than linearisation data treatments in the complexing studies aimed to define behaviour of ligands in natural waters.

Those studies are based on experimental data obtained from sample titrations, i.e. metal additions, coupled with ASV, CSV or ISE measurements [26,28,[34][35][36]. The number of standard additions, their repartition and covered concentration range, which characterise these titrations, are generally quite arbitrary. However, for the analysis of multi-ligand systems, these characteristics are certainly important to correctly determine the complexing parameters, especially since they are dependent on the analytical technique used (i.e. its detection limit, analytical window, etc.).

In the present work, we analyse the influence of metal standard addition repartition (linear, decade and logarithmic) on the results of simulated experimental data fitting. The same simulated experimental data points will be treated with the methods of Chau-Buffle [5,9], Scatchard [8], Ružić-van den Berg [10], non-linear fitting [32] and PROSECE [37]. The objective of the present work is to optimise experimental conditions to obtain the most reliable results. The work has determined (i) the uncertainty of metal complexing capacity parameters calculated by each method, and (ii) which method is the most rigorous for each simulated L iT , K i example.

Simulation of "experimental" data points

Simulation of experimental data points, i.e. their numerical generation, is useful for the comparative analysis that will be provided in this work, where we focus on a simple case of a 1:1 metal-ligand system and analyse the cases of one-and two-ligand complexes with a single metal. To calculate each point according to the design of experimental standard additions, the program MINEQL [START_REF] Westall | MINEQL: a computer program for the calculation of chemical equilibrium composition of aqueous systems[END_REF] is used. It is not a trivial calculation, as given by Pižeta and Branica [32] to find the exact-metal free concentration from the metal added, when the number of ligands is larger than 1.

The concentration ranges and the number of additions of titrated metal were determined as a function of the analysed ligands. The range limits were based on the analytical window of common electrochemical techniques (ASV or CSV) used when measuring metal complexing properties of natural samples. It has been stated in the literature that it was difficult to correctly determine more than four dif-ferent ligand groups [18,20,34]. In practice, the concentration range very rarely covers more than four decades (e.g. pM T = -log(M T ) from 10 to 6, where M T is total metal concentration). Application of many standard additions is facilitated by automatic burettes controlled by a computer, however, this is limited in practice by the total duration of the titration. For example, one titration of Suwannee River fulvic acid with cadmium and lead (30 additions in a logarithmic mode) using differential pulse (DP) ASV (Metrohm stand 663 + Autolab + 4 Metrohm burettes) required more than 48 h for complete characterisation of the complexing properties [START_REF] Garnier | [END_REF]. According to this knowledge we estimated as follows: for one ligand the concentration range of metal additions can be fixed to two decades (1-100 nM) defined by 20 points considered to be a minimum for complexing studies, and for four ligands it can be fixed to four decades (1-10,000 nM) by 35 points, which are a priori realistic and sufficient, compared to the experiments performed with natural samples [15,17,19,36]. By interpolation, a simple calculation gives 2.66 decades for two ligands, i.e. 1-457 nM, defined by 25 standard addition points.

To analyse the effect of titration mode on simulation, three different strategies of non-arbitrary standard additions of metal ions were simulated: linear, decade and logarithmic additions, bearing in mind the possibility of utilising automated burettes in a real experiment.

These three modes consist, respectively, from i = 1 to n additions:

• constant increments of M T , i.e. M T,i = M T,n /n, • constant increments of M T in one decade of concentration, i.e. for 1-10 nM: M T,i = i×10 -9 , for 10-100 nM: M T,i = (i -10) × 10 -8 , etc., • constant increments of pM T , i.e. M T,i = 10 -(pM T,0 -(pM T,0 -pM T,n )×(i/n)) , where n is the total number of points in the titration, M T,0 , M T,n , M T,i are the initial, final and ith total metal concentrations.

In order to obtain realistic data, taking into account previous results of Pižeta and Branica [32] about the influence of noise added to the data (0-10%) and considering an increase in the sensitivity of electrochemical instruments, a random noise of 2% was added to the points after the simulation (MINEQL) and before fitting. The set of applied random noise values varies between two simulated experiments with different L iT , K i values, but remains unchanged for one experiment treated by the three titration modes.

The position of the points selected for simulation in log K × log L T space, compared to Pižeta and Branica [32], is presented in Fig. 1.

Methods of data treatment and fitting to the model

The methods of data treatment for metal complexing capacity and stability constant determination comprise three linearisation methods: (a) linear regression proposed by Chau [5] in combination with linearisation proposed by Buffle [9], (b) linearisation by a mathematical transformation proposed by Scatchard [8], and (c) Ružić [10] and van den Berg [11] followed by linear regression; as well as two optimisation methods: (i) non-linear fitting of the data transformed according to Ružić-van den Berg method when they are not linear after the transformation, proposed by Pižeta and Branica [32], and (ii) a matrix approach by the program named PROSECE (summarized in Table 1). M T consists of the metal initially present in the solution (M 0 ) and the metal added (M a ). M f means free metal, not inertly bound to the complexing ligand, which is a measured value in voltammetry. Theoretically, the methods based on mathematical transformation of data (linearisation in this case) are more sensitive to the noise of M f than the methods using non-modified data. Indeed, transformations induce noise propagation which influences the two axes y,

x = f(M f , M T ).
PROSECE (Programme d'Optimisation et de SpEciation Chimique dans l'Environnement), was developed in RCMO-PROTEE laboratory (with the help of P. Seppecher (ANLA laboratory, University of Toulon and Var)), for complexing studies using a discrete model. It has been written on the basis of the numerical calculation program Octave (freeware). This new software is available on request. It uses a speciation subroutine coupled to an optimisation one. The speciation calculation is based on mass balance resolution by the Newton-Raphson optimisation method, after Morel's table definition representing all complexing equations between the principal components (M, L i , H, etc.) and the component species and compounds (M f , L if , ML i , etc.) in the defined chemical system [START_REF] Morel | Principles of Aquatic Chemistry[END_REF]. After the definition of the initial values for complexing parameters (L iT , K i ) which must be determined, the speciation modulus calculates the concentrations of component species and compounds. Comparison with experimental (or simulated) data gives an error value. The optimisation process, by means of a modified simplex [START_REF] Massart | Handbook of Chemometrics and Qualimetrics: Part A[END_REF], adjusts these selected values of complexing parameters, giving a new set of parameters for the speciation part, which continues until minimum error value has been reached. For instance, PROSECE allows the modelling of natural organic matter complexing and competing properties versus proton and trace metals by optimisation of a set of complexing parameters (total ligand concentrations, acidity constants, stability constants versus trace metals) [37].

One-ligand systems

From previous knowledge and experience of the detection and determination of metal complexing capacity with one ligand, a two decades range of metal standard additions is necessary. As Ružić [10] has shown that K × L T > 1 is necessary to obtain reliable results, the analyses and comparisons were made for the case K × L T = 50 in the pairs as follows: L T = 10, 50 and 90 nmol l -1 and K = 5 × 10 9 , Table 1 Overview of the data treatment methods for metal complexing capacity determination

Graphical presentation

Mathematical treatment and explanation One-ligand system: linearisations Linear regression of last points of the curve M f = aM T + b, L T is determined by extrapolation:

L T = -b/a, linearisation of experimental points according to L T /(M T -M f ) = (c/M f ) + d, where L T is determined by the method of Chau. K = 1/c (d = 1)
Mathematical transformation of experimental data points according to

(M T -M f )/M f = f(M T -M f ), linear regression: (M T -M f )/M f = a(M T -M f ) + b, K = -a, L T = -b/a
Mathematical transformation of experimental data points according to

M f /(M T -M f ) = f(M f ), linear regression: M f /(M T -M f ) = aM f + b, K = a/b, L T = 1/a
Two-ligand system: linearisations Mathematical transformation of experimental data points according to

(M T -M f )/M f = f(M T -M f )
, linear regression of the first data points and the last ones:

Y 1 = a(M T -M f ) + b and Y 2 = c(M T -M f ) + d, calculation of L 1T , L 2T , K 1 , K 2 = f(a, b, c, d) Table 1 (Continued )

Graphical presentation Mathematical treatment and explanation

Mathematical transformation of experimental data points according to

M f /M T -M f = f(M f ), linear regression: Y 1 = (M f /(M T -M f )) M f →∞ = aM f + b, Y 2 = 1/(aM f + b -(M f /(M T -M f )) M f →∞ ) = cM f + d, calculation of L 1T , L 2T , K 1 , K 2 = f(a, b, c, d)
One-and multi-ligand systems: fitting After Ružić-van den Berg linearisation, for more than one-ligand (N > 1), non-linear fitting based on Levenberg-Marquardt method:

M f /(M T -M f ) = 1/ N i=1 (L iT /M f + (1/K i ))
Speciation calculus based on mass balance table resolution. Comparing the "experimental" and the calculated concentrations, an optimisation process, using a modified simplex, adjusts the complexing parameters until reaching an optimum error 1 × 10 9 and 5.56 × 10 8 l mol -1 , respectively. In Fig. 1, it is presented that the points are in the probable and convenient areas of possible experimental situations. Simulated data points, representing each pair of L T and K for different types of standard additions in linear, decade and logarithmic modes are presented in Fig. 2.

All the results for a one-ligand system are summarised in Tables 2 and3. The values obtained by data treatment are given together with the preset values (K th , L th ). When it was possible to calculate these, the values of confidence interval (e 1 ) of the calculated parameters are given in % ( K% and L T %). The difference between the preset and the calculated values is called e 2 . Correlation of log(e 1 ) versus log(e 2 ) is summarized in Table 3. It shows an important correlation for linearisation methods, which could be used in predicting the efficiency of parameter determination.

Linear additions

The results obtained for linear mode have drawn clear distinction between different techniques of data treatment. Thus, the errors for K and L T fluctuate from 0.1 to 122.9%, and from 0 to 54.5%, respectively, depending on the experiment and on data treatment. PROSECE is the only method which allows correct K and L T determinations, whichever the L T value in the studied analytical window (errors between 0.6 and 3.9% for K, and between 0.0 and 0.6% for L T ).

As long as the L T value is large enough, correct results are obtained for Ružić-van den Berg and Scatchard lineari- Chau-Buffle linearisation gives reliable results only when L T is centred on the analytical window (shape of experiment 2 in Fig. 2). Indeed, a small preset value for L T leads to an incorrectly calculated L T value since only the last few titration points were used for Chau linearisation. Further, when those titration points are too far from the cure (i.e. complexing capacity equivalent value) they are sensitive to the noise (shape of experiment 1 in Fig. 2). A large value for L T also leads to a large error for L T , because M T,n is small compared to the theoretical L T , which is leads to the situation that the last titration points are in the curved portion, i.e. narrowly reaching complexing capacity equivalent value (shape of experiment 3 in Fig. 2). These two different phenomena increase the error for K determination by Buffle linearisation using Chau L T value.

Generally, a linear titration seems to be adapted to L T and K determinations only for high ligand concentrations when there are enough data before the curved portion. In that case, Ružić-van den Berg and Scatchard linearisations give reliable results. However, data treatment by PROSECE uses all data points without mathematical transformation and achieves correct results even for small L T values.

Decade and logarithmic additions

In general, the results obtained by decadic and logarithmic additions are improved compared to the linear additions results, but they do not allow these two non-linear titration modes to be distinguished.

Scatchard linearisation and PROSECE modelling allow the determination of K and L T with very low error values (<1.7%), regardless the L T position on the metal range scanned. Thus, these titration modes increase the sensitivity for complexing parameters optimisation, especially for small values of L T .

In contrast, they do not increase the Ružić-van den Berg ability to correctly determine the stability constant, which is always negative when the L T value is low.

The errors for K and L T with Chau-Buffle linearisation are slightly improved, particularly in the case of low L T values.

For the case of a 1:1 metal-ligand system, when carried out at constant K×L T values, these experiments have shown that a fitting process by using non-modified data, such as PROSECE, is the most efficient technique, despite the apparent simplicity of the system. Indeed, the results obtained are close to the predetermined values, regardless of the titration mode or the L T value used.

By better defining low metal concentration ranges, logarithmic and decade additions greatly increase the sensitivity of Scatchard linearisation for the determination of the complexing parameters whichever the ligand concentration. In contrast, Ružić-van den Berg linearisation, dependent on high metal concentrations (due to the mathematical transformation used), seems to be less adapted to the stability constant determination when the ligand concentration is small, which happens frequently in natural systems. The errors obtained are not minimised by the use of non-linear titrations.

The Chau-Buffle method, using two successive linearisations, produces less reliable results than those obtained by the other methods.

Two-ligand systems

In contrast to one-ligand case, where only three simulated experiments have been carried out with the preset K × L T values, a factorial matrix [START_REF] Massart | Handbook of Chemometrics and Qualimetrics: Part A[END_REF] has been used for two-ligand systems to evaluate the number and the type of experiments to be performed. Considering two ligands defined by four parameters (K i and L iT ), 2 4 = 16 experiments have been simulated by using high (+) and low (-) values for each parameter. Previous results, showing the influence of the values of the complexing parameters and of random noise added to the retrieving of the simulated data [10,14,18,32,33], as well as the use of different analytical techniques and different models to determine those parameters in natural samples [15][16][17][18][19][20][21][22][23][24][25][26][27][28][29][30], led us to choose the K i and L iT values for this study. Thus, (+) and (-) values for K 1 , L 1T , K 2 , L 2T are, respectively, 10 10 and 10 8 , 30 and 10 nM, 10 7 and 10 6 , and 300 and 100 nM (i.e. K 1 × L 1T : 1-300, K 2 × L 2T : 0.1-3, K 1 /K 2 : 10-10 4 ). Table 4 shows the values used for these 16 simulated experiments, using Yates matrix design [START_REF] Massart | Handbook of Chemometrics and Qualimetrics: Part A[END_REF]. Three of these simulated experiments (n • 3, 5 and 7) have been plotted using M f = f(M T ) and pM f = f(pM T ) scales, in linear, logarithmic and decade titration modes (Fig. 3). Experiments 3 and 5 could not be distinguished in linear mode, whatever the considered scale, but were significantly different in logarithmic or decade titration modes using the pM f = f(pM T ) scale.

To each simulated experiment, random noise of 2% was added to the calculated points of free metal concentrations (M f ). In order to be better defined and not biased by the knowledge of the correct answers, optimisation methods (non-linear fitting and PROSECE) always employ the same initial values for the calculated parameters (K 1 = 10 9 , L 1T = 20 nM, K 2 = 3.16 × 10 6 , L 2T = 200 nM), which are not necessarily optimal, i.e. in the case of real experimental data treatment, one would change the initial parameters in order to optimise the results, especially if negative results occur. In order correctly to compare the efficiency of titration modes (linear, logarithmic and decade) and of the methods for data treatment (Ružić and Scatchard linearisations, non-linear fitting and PROSECE), errors for four different parameters have been calculated. They are based on the differences between the theoretical value of complexing parameter (X th ) and the calculated value (X calc ), according to the following equation:

Linear additions

X error = X th -X calc 2
X th 2 × 100.

Theoretically, using Yates method, it would be possible to quantify the effects of these different factors (i.e. K i and L iT ) on the responses (i.e. errors values) with the use of equations which would allow to determine the influence of K i and L iT values on the data treatment. However, this is based on the assumption that response varies linearly when a parameter changes from (+) to (-) value. In our case, it does not seem realistic, and the calculations have shown that it is impossible to link modelling error with the preset values of K i and L iT .

To estimate possible correlations between these different parameters, the obtained errors have been presented as:

K 1 error = f(L 1T error ), K 2 error = f(L 2T error ), L 2T error = f(L 1T error ), K 2 error = f(L 1T error ), K 1 error = f(L 2T error )
and K 2 error = f(K 1 error ). Thus, the error for K 1 seems to be well correlated to the L 1T error, which emphasizes the reciprocal influence of strong ligand binding parameters. The same phenomenon, although less pronounced, is obtained between K 2 and L 2T errors. In contrast, K 2 versus L 1T and K 1 versus L 2T errors have shown inverse correlations, i.e. wrong determination of one parameter is often linked to the correct determination of the other one, and vice versa.

However, except for those tendencies, relatively few correlations can be observed, which shows the complex dependence between all these complexing parameters.

Linear additions

Undeniably, PROSECE fitting gives the best results for K i and L iT determinations, independently of the preset values. Only a few of the simulated experiments led to significant errors for L 1T and L 2T . These corresponded to the experiments where the preset parameter values were low.

In contrast, Ružić linearisation with its mode of addition does not allow a correct determination of K 1 , to whichever of the preset values applied. Indeed, for 10 of them, the obtained values are negative. Retrieving of ligands concentrations is correct, almost independent of the conditions (variations in errors are small for all experiments) but less than for PROSECE. Determinations of K 2 are always better than those of K 1 ; they are more accurate the K 2 × L 2T becomes larger.

For the determination of K 1 , Scatchard linearisation seems to be markedly more convenient than Ružić method. However, the results obtained by this linearisation do not reach the reliability of those obtained by PROSECE. In contrast, in retrieving L 2T with negative values for two experiments, it is less accurate than Ružić method. This shows the importance and the influence of mathematical transformation on the final result, i.e. the values of complexing parameters.

Despite data linearisation when using Ružić method, non-linear fitting seems to be less sensitive with regard to K 1 calculation. However, with five negative values for K 1 and two negative values for L 1T , the obtained results differ from the preset values. Also, K 2 and L 2T errors are often larger than those obtained by Ružić and Scatchard linearisations.

Linear additions do not seem to be adapted for determination of complexing properties of strong ligands at low concentrations, which are frequently found in natural waters, except if PROSECE is used. With PROSECE the errors are 30% on average compared to other methods which gave errors of 50-500%. This can be probably explained by an insufficient definition of the concentration range between 1 and 50 nM (only three points).

However, as regards the complexing properties of the weaker ligands of higher concentrations, despite large number of data defining them (20 data points from 100 to 457 nM), linearisations and non-linear fitting do not allow accurate determination of K 2 , possibly due to too low L 2T × K 2 values. Only L 2T is correctly determined (20% error on average).

Logarithmic additions

Results obtained for the determination of K 1 and L 1T are clearly improved with the use of logarithmic additions compared to linear ones, whatever the data treatment used.

This titration mode allows a complete elimination of all the aberrant results occurred in Ružić linearisation, i.e. negative K 1 values. Except for one simulated experiment, the same happened to the results from non-linear fitting. Thus, mean errors for L 1T , for these two methods of data treatment, change respectively from 43.9% (linear) to 11.4% (logarithmic) for Ružić linearisation and from 512.4 to 15.6% for non-linear fitting. So, mean errors for K 1 decrease, respectively, from 1457 to 8.5% for Ružić linearisation and from 96.6 to 52.6% for non-linear fitting.

Except for one simulated experiment with an L 1T error of 220%, PROSECE allows determination of complexing parameters of the first ligand with mean errors <9% (50% of the experiments with an error of 3%). It seems to be the most suitable method for the determination of complexing properties of low concentrations of strong ligands.

Although markedly improved, results obtained by Scatchard linearisation are not in the same range (mean errors: 20%).

Concerning L 2T determination, linearisations and nonlinear fitting give slightly more correct results than PROS-ECE, which is influenced by the values of K 1 /K 2 and those of L 1T /L 2T . For the same parameter, Scatchard linearisation benefits the most when changing from linear to logarithmic mode.

The K 2 parameter remains better determined by PROS-ECE than by linearisation data treatments, which give very large errors. Non-linear fitting and linearisation data treatments, which are very sensitive to L 2T and K 2 values compared to the L 1T and K 1 ones, lead to significant errors for a large number of experiments.

Generally, non-linear fitting and linearisation data treatments will correctly determine the concentration of weak ligands, but only PROSECE can determine the stability constant of these ligands, despite the low values of K 2 × L 2T (between 0.1 and 3).

Logarithmic additions, improving markedly the titration accuracy at low concentrations (16 points from 0 to 50 nM, compared to 3 points in the linear mode), amplify the resolution of the complexing properties of low concentration of strong ligands regardless of the data treatment used. This phenomenon is obvious for K 1 determination by Ružić linearisation and non-linear fitting.

Despite a decreasing number of points defining high concentrations (7 points from 100 to 457 nM compared to 20 points for linear additions), the results of Ružić linearisation and PROSECE determinations of the complexing properties of weaker ligands are similar to those obtained in the linear mode. They are even slightly improved as regards K 2 treated by non-linear fitting and L 2T calculated by Scatchard.

Decade additions

The results obtained with this titration mode are in general very close to those calculated from the logarithmic ones. So, these additions largely improve the determination of complexing parameters for low concentrations of strong ligands compared to linear additions.

The resolution obtained by Ružić linearisation seems to be slightly better than that by logarithmic additions. This data treatment appears to be suitable for the calculation of K 1 and L 1T but only when the precision in the analytical window is sufficient at low total metal concentrations.

Although less correct, Scatchard linearisation leads to errors for L 1T and K 1 < 15%.

With one simulated experiment leading to a negative K 1 , and mean errors of 25%, non-linear fitting seems to be the least adapted of the four data treatments.

Despite one simulated experiment with large error for K 1 with PROSECE, errors are slightly improved compared to logarithmic additions (mean error 5%, except the maximal error for K 1 ). Similar to other addition modes, this data treatment seems to be the most adapted to determine complexing of strong ligands. In contrast to L 1T and K 1 , determination of L 2T and K 2 is not improved by the use of this titration mode compared to the logarithmic mode. Linearisations give correct value for L 2T but poor results for K 2 . Non-linear fitting improves a little the determination of the latter. Similar to logarithmic additions, PROSECE seems not to be the best adapted to define L 2T correctly but it is the only one that allows an almost correct calculation of K 2 .

Similar to the logarithmic mode, this titration mode, by largely increasing the precision of the lowest concentrations (16 points compared to 3 points in the linear case), improves the resolution of complexing parameters for strong ligands.

In contrast, the number of points in a higher concentration range being smaller than for other titration modes (only five points between 100 and 457 nM), leads to poorer determination of weaker ligands.

Influence of random noise distribution

In order to approach realistic [M f ] values, random noise between -2 and 2% was added to the values simulated by MINEQL. However, it seems necessary to comprehend the effect of noise distribution on the complexing parameters determination. For this purpose, the data from two previous simulated experiments (no. 5 in linear mode and no. 7 in logarithmic mode, cf. Table 4) have been treated by all methods, after the addition of 10 different sets of random noise. The first simulated experiment corresponds to the one for which it was difficult to retrieve the preset values in all data treatment methods, especially with linearisations and non-linear fitting that have led to negative K 1 , L 1T and/or L 2T values. In contrast, the second experiment gave correct results whichever data treatment method was used.

Theoretically, linearisations should be more sensitive to these random noise distributions than fitting (especially for the first experiment), because, in contrast to fitting, they take into account only a part of the data points which increase the noise effect. Moreover, mathematical transformations by data linearisation affect X and Y values (X, Y = f([M f ])), while the use of non-transformed data as in the case of PROSECE affects only Y values (cf. Table 1). All obtained results are summarised in Table 5.

As expected, the techniques based on data transformation (followed by linearisation or fitting) appeared to be more sensitive to data noise (i.e. mean values were far from the preset ones, with large standard deviations, and a great number of negative values) than PROSECE, which is working with non-modified data. Thus, this data treatment can lead to aberrant values (e.g. negative values) of the parameters of the first and/or the second ligand, but depending on random noise distribution, even more so as K 2 × L 2T factor or L 1 become low, which is a common situation in natural conditions. When K 1 × L 1T and K 2 × L 2T are large, i.e. greatly larger than 1, the four data treatment methods give quite similar results. Fig. 4. "Contour-graphs" of "K 1 error = f((L 1T × K 1 ), (L 2T × K 2 ))" obtained with Ružić linearisation in linear (A) and logarithmic (B) additions ((+) calculated K 1 error , () calculated but negative K 1 , (ᮀ) calculated and higher than 1000).

Possibility of error prediction

Taking into account the results obtained in these 16 experiments, with different K 1 × L 1T and K 2 × L 2T values, the question arose whether it would be possible to predict errors of the complexing parameters for simulated experiments with other preset values and to estimate the accuracy of the results obtained by treating the data of a real experiment obtained from natural samples analysis (Fig. 3).

Considering the K 1 × L 1T and K 2 × L 2T values of 16 experiments, "contour-graphs" showing parameter errors = f((L 1T × K 1 ), (L 2T × K 2 )) have been plotted in a logarithmic scale for each couple of titration mode-data treatment method, using Surfer 6 software. From a large number of graphs obtained (48), only two are presented in Fig. 4. All figures could be found at: http://mounier.univtln.fr/ACAgarnier2003.html.

To test the possible prediction capability of these graphs, five new experiments, with different K i × L iT values, were simulated in linear and logarithmic modes. To estimate er- L 1T (nM) 4 and the parameters of these new simulated experiments are summarised in Table 6. These six experiments were treated using all the defined techniques to obtain calculated error values. For each experiment, error estimation consists of reading-out, on a "contour-graph", the error value at the K i × L iT theoretical point, and comparing it with the value actually obtained.

K 1 L 2T (nM) K 2 L 1T × K 1 L 2T × K 2 K 1 error linear (%) K 1 error logarithmic (%) 1
Theoretically, three different cases are expected to be obtained depending on whether the estimated values are:

• Close to the calculated ones: graphs are adapted to predict accuracy of a given titration mode-data treatment couple, to determine complexing parameters. • Higher than the calculated ones: prediction is pessimistic, which is not a disadvantage since graphs allow estimation of upper limit values.

• Inferior to the calculated ones: the proposed plots are not adapted to estimate the complexing parameters errors for a non-simulated experiment.

By plotting the estimated errors versus the calculated ones, for linearisation techniques, correlation is obtained. It confirms that "contour-graphs" allow the correct estimation of error (Fig. 4B), except when the obtained parameter values are incorrect (such as a negative K 1 obtained by the Ružić method in linear mode (Fig. 4A)).

In contrast, predictions for the fitting techniques are less correct. In most cases prediction values overestimate the real values. Despite the important gap obtained, predicted errors can even be used as upper limits.

Regarding the variability of errors for different parameters values, with the same K i × L iT values (i.e. the difference between experiments 3 and 6), in all cases, standard deviations of calculated errors between these two experiments were >50%. Moreover, differences between the predicted values and the mean calculated values are, in general, smaller than those standard deviations. This means that, even if the variations of errors in parameters determinations are large, for different parameters with combinations with the same K i × L iT values, the errors predicted with "contour-graphs" are in the variation domain of real errors, and are therefore usable for prediction.

Conclusions

The provided analyses of the influence of the titration type and of data treatment method on the accuracy of metal complexing parameters determination for one-and two-ligand systems based on the simulated data corresponding to those obtained by anodic stripping voltammetry, suggest the following conclusions:

Even with one ligand system, logarithmic and decadic titrations are much better than the linear ones.

Linearisation methods are in many cases inferior to the optimisation methods, especially in a linear titration mode. Random noise has a significant influence on the results of linearisation methods (especially in the case of two-ligand systems), while only a part of the data points participates in the determination of complexing parameters. On the other hand, for linearisation methods, in the case of one-ligand systems, a high correlation between the confidence interval of the calculated parameters and the difference between the preset and the calculated values is found. This correlation is proposed to be used as an estimate for the quality of the results in real experiments.

PROSECE is by far superior to other methods (linearisations and non-linear fitting) in most cases due to its flexible and powerful mathematical background. It is a result of computer programming development and is highly recommended as a tool for data treatment. Using non-transformed data, PROSECE minimizes the effect of data noise prop-agation, which seems to be favourable for complexing parameters determination compared to other data treatment methods defined.

Construction of "contour-graphs" enables certain error prediction of the calculated complexing parameters and is proposed as an orientation and valorisation tool in real samples analyses.

Fig. 1 .

 1 Fig. 1. Repartition of ligands complexing properties used in this study ((᭹) one-ligand, (᭺) first of two-ligand, (᭛) second of two-ligand) and in Pižeta and Branica previous study ((᭡) one-ligand, (×) first of two-ligand, (ᮀ) second of two-ligand).

Fig. 2 .

 2 Fig. 2. No-noised simulated additions for one ligand model experiments (see the corresponding preset parameters values in Table 2 marked by (a)).

Fig. 3 .

 3 Fig. 3. No-noised simulated additions for two ligands simulated experiments nos. 3, 5 and 7.

Fig. 4 .

 4 (×) No.

Table 2

 2 Results obtained for one ligand system

		Linear			Logarithmic			Decade
		Exp. 1	Exp. 2	Exp. 3 a	Exp. 1 a	Exp. 2	Exp. 3	Exp. 1	Exp. 2 a	Exp. 3
	K th	5.00E+09	1.00E+09	5.56E+08	5.00E+09	1.00E+09	5.56E+08	5.00E+09	1.00E+09	5.56E+08
	L Tth	1.00E-08	5.00E-08	9.00E-08	1.00E-08	5.00E-08	9.00E-08	1.00E-08	5.00E-08	9.00E-08
	K	3.95E+09	1.05E+09	7.25E+08	4.50E+09	1.02E+09	7.62E+08	4.01E+09	1.06E+09	6.94E+08
	Chau-Buffle linearisations
	K (%)	110.7	0.5	0.6	1.3	0.6	0.7	1.3	0.5	0.5
	L T	1.54E-08	4.75E-08	6.88E-08	1.10E-08	4.87E-08	6.54E-08	1.23E-08	4.70E-08	7.19E-08
	L T (%)	121.4	18.7	23.5	33.4	ND	ND	0.0	0.0	ND
	K	-1.15E+09	1.03E+09	5.58E+08	-2.54E+09	1.02E+09	5.56E+08	-3.55E+09	1.01E+09	5.56E+08
	Ružić-van den Berg linearisation
	K (%)	-466.8	11.9	1.0	-427.5	5.8	0.6	-779.4	7.7	0.7
	L T	9.48E-09	4.96E-08	8.97E-08	9.25E-09	4.95E-08	8.96E-08	9.45E-09	4.96E-08	8.97E-08
	L T %	7.7	0.5	0.2	5.4	0.4	0.2	5.6	0.4	0.2
	K	9.52E+08	1.00E+09	5.55E+08	5.00E+09	9.97E+08	5.57E+08	4.98E+09	9.99E+08	5.56E+08
	Scatchard linearisation
	K (%)	67.2	0.0	0.0	0.2	0.0	0.0	0.3	0.0	0.0
	L T	1.04E-08	4.97E-08	8.99E-08	9.90E-09	4.98E-08	8.96E-08	9.90E-09	4.98E-08	8.97E-08
	L T (%)	162.4	1.5	1.3	7.6	1.3	1.5	9.7	1.2	1.3
	K	4.81E+09	1.01E+09	5.62E+08	4.92E+09	9.90E+08	5.59E+08	4.94E+09	1.00E+09	5.59E+08
	PROSECE
	K (%)	ND	ND	ND	ND	ND	ND	ND	ND	ND
	L T	1.00E-08	4.97E-08	8.95E-08	1.00E-08	4.98E-08	8.96E-08	1.00E-08	4.97E-08	8.96E-08
	L T (%)	ND	ND	ND	ND	ND	ND	ND	ND	ND

a Simulated data shown in Fig.

2

.

Table 3

 3 Correlation coefficient (R) between log(e 1 ) and log(e 2 ) values

					sations (errors <3% for K and L T ). In contrast, the first ex-
					periment, representing low concentration ligands, leads to
					negative K values for Ružić-van den Berg and high error for
					Scatchard, showing unsuitability of the linear mode treated
					by linearisations for the lowest L T values.
		Chau-Buffle	Ružić-van	Scatchard	PROSECE
		linearisations	den Berg	linearisations
			linearisations
	K 1	0.231	0.993	0.864	ND
	L 1T	0.560	0.980	0.919	ND

Table 4 The

 4 Case of two ligands, preset complexing parameters values for simulated experiments

		+	-	No.
				1	2	3	4	5	6	7	8	9	1 0	1 1	1 2	1 3	1 4	1 5	1 6
	L 1T	30 nM	10 nM	+	+	-	+	-	-	+	+	-	-	+	+	-	+	-	-
	K 1	10 10	10 8	-	-	-	-	+	-	+	+	-	+	+	-	-	+	+	+
	L 2T	300 nM	100 nM	-	+	-	+	-	-	+	-	+	+	-	-	+	+	+	-
	K 2	10 7	10 6	-	-	-	+	-	+	+	-	-	+	+	+	+	-	-	+

Table 5

 5 The influence of random noise distribution on the complexing parameters determination for different data treatment methods (mean values, standard deviation and number of negative values calculated for 10 different sets of noise-added data in simulated experiments nos. 5 and 7)

		Ružić linearisation		Scatchard linearisation		PROSECE			Non-linear fitting
		Mean	S.D.	N < 0 Mean	S.D.	N < 0 Mean	S.D.	N < 0 Mean	S.D.	N < 0
		value	(%)		value	(%)		value	(%)		value	(%)
	Exp. 5: linear theory
	L 1T 1.00E-08 5.65E-09	27.6	0	1.44E-08	3.9	0	1.06E-08	5.3	0	5.66E-09 100.6	4
	K 1	1.00E+10 -2.94E+08 109.9	10	3.08E+08	18.5	0	5.01E+09 138.1	0	5.83E+08 215.9	5
	L 2T 1.00E-07 8.71E-08	51.2	0	-1.51E-07 294.4	7	2.43E-07	20.0	0	5.05E-07 285.0	1
	K 2	1.00E+06 2.32E+06	45.1	0	-2.33E+05 327.1	7	3.71E+05	31.5	0	4.25E+06	93.4	1
	Exp. 7: logarithmic theory
	L 1T 3.00E-08 2.99E-08	2.4	0	3.06E-08	2.2	0	3.00E-08	0.4	0	3.00E-08	3.0	0
	K 1	1.00E+10 1.00E+10	2.9	0	9.77E+09	2.6	0	9.99E+09	0.7	0	1.01E+10	7.6	0
	L 2T 3.00E-07 2.83E-07	3.1	0	2.80E-07	3.1	0	2.95E-07	2.4	0	2.97E-07	2.6	0
	K 2	1.00E+07 1.22E+07	6.6	0	1.38E+07	6.8	0	1.03E+07	3.8	0	1.02E+07	6.0	0

Table 6

 6 Parameters values of the additional simulated experiments for error prediction test and the errors obtained for K 1 , for the data treatment by Ružić linearisation for linear and logarithmic titration modes (errors values in italic correspond to the negative values of K 1 ), viz.

Acknowledgements

The authors wish to thank P. Seppecher (ANLA laboratory, University of Toulon and Var) for his substantial contribution to the development of PROSECE, as well as the PhD school of the University of Toulon and Var for financial support of collaboration between France and Croatia.


C Garnier

S Mounier

K Küpker

Y Lucas

J Benaïm

MODELISATION DU COMPORTEMENT DE LA MON PAR UN « CHIMIO-TYPE »

Keywords: MON, métaux traces, proton, DPASV, PROSECE NOM, trace metals, proton, DPASV, PROSECE

L'analyse du comportement d'acides fulviques et humiques standards (Suwannee River et Laurentian River), vis-à- vis des métaux traces (Cd, Pb) et du proton, a été réalisée à partir de différentes titrations potentiométriques et voltammétriques, par additions de solutions standards en mode logarithmique. Les données expérimentales ainsi obtenues ont été traitées à l'aide du logiciel PROSECE, en définissant un modèle basé sur un « chimio-type » constitué de propriétés discrètes analysables des MON étudiées. Ce « chimio-type », constitué d'un set de trois types de « quasi-particules » (Sposito, 1981), a permis de modéliser les sites acides (6 « quasi-particules » de type I), les sites complexants échangeables (4 de type II) ainsi que les sites complexants spécifiques (2 de type III) présents sur la MON. Les paramètres de complexation optimisés (i.e. concentrations, constantes de complexation et d'acidité de ces « quasi-particules ») permettent de prévoir le comportement de la MON vis-à-vis des métaux traces, et les phénomènes de compétition métal/métal et métal/proton. Ces paramètres, quasithermodynamiques, sont facilement intégrables dans un modèle de transport de contaminants tels que le modèle SiAM3D-MOCO (IFREMER).

INTRODUCTION

Dans les milieux naturels, tel un estuaire, les composés constituants la matière organique naturelle (MON) ont un rôle prépondérant quant à la complexation des métaux traces [START_REF] Buffle | Complexation reactions in aquatic systems[END_REF]. En effet, la spéciation chimique des métaux traces dans les eaux naturelles est fortement liée à la composition chimique de l'environnement (i.e. concentrations en éléments majeurs, autres métaux, proton, MON, ligands organiques…). De ce fait, la MON contrôle en partie le transport de ces polluants métalliques ainsi que leur réactivité vis-à-vis du milieu, notamment des micro-organismes (phytoplancton, zooplancton, macro-algues, …), en conditionnant la toxicité et/ou la biodisponibilité des métaux [START_REF] Vasconcelos | Antagonistic interactions of Pb and Cd on Cu uptake, growth inhibition and chelator release in the marine algae Emiliania huxleyi[END_REF]. L'étude de la MON, en vue d'une modélisation de ses propriétés complexantes et en conséquence d'une prédiction de son comportement dans un milieu donné, est donc un enjeu majeur dans la compréhension de l'impact des métaux traces et de la MON elle-même dans l'environnement. Ainsi, la caractérisation de la MON par des paramètres de complexation, à l'aide de modèles, a fait et fait encore l'objet de nombreux travaux [START_REF] Town R | Dispelling the myths: Is the existence of L1 and L2 ligands necessary to explain metal ion speciation in natural waters?[END_REF]. Compte tenu de la diversité des techniques d'analyses utilisées (ASV, CSV, ISE, fluorescence…) et des modèles appliqués (continus ou discrets), les données de complexation obtenues sont difficilement comparables et utilisables dans un modèle hydrogéochimique (SiAM3D-MOCO, IFREMER, [START_REF] Gonzalez | Modeling of cadmium speciationand dynamics in the Seine estuary (France)[END_REF]. Le travail ici présenté a pour objectif de caractériser des MO extraites, par un « chimio-type » constitué d'un ensemble de paramètres de complexation, vis-à-vis des métaux traces (Cd, Pb) et du proton, directement utilisables dans un modèle de transport des contaminants. Afin de déterminer une distribution de ligands sur la MON la plus grande possible (i.e. de 0.1nM à 10µM de [M] T ), deux techniques d'analyses ont été utilisées (voltammétrie, potentiométrie). Pour augmenter la précision sur la détermination des propriétés de complexation et de compétition, une technique de titration par incréments constants en logarithme de la concentration totale en métal ([M] T ), a été développée. De même, un programme d'optimisation et de spéciation chimique (PROSECE) a été mis au point afin de traiter les données expérimentales, ceci dans l'objectif de prédire le comportement de la MON grâce à un set de paramètres quasi-thermodynamiques.

MATERIELS ET METHODES

Matières organiques utilisées

Les solutions de MO extraites utilisées dans cette étude ont été réalisées à partir d'Acides Fulviques (AF) et Humiques (AH) Suwannee River (SR) et Laurentian River (LR). La force ionique des solutions utilisées a été ajustée à 0.10M par NaNO 3 (J.T. Baker). Les titrations électrochimiques ont été réalisées pour des teneurs en Carbone Organique Dissous (COD) proches de 10mg/L, et en utilisant des tampons acétique et borax (théoriquement peu complexants vis-à-vis des métaux) afin d'ajuster le pH, respectivement à 4.6 et 7.8. Les solutions sont stockées à 4°C dans le noir. Les titrations potentiométriques ont quant à elles été faites avec des teneurs en COD allant de 100 à 200ppm de COD, concentrations équivalentes à celles utilisées dans des études similaires (Ritchie & Perdue, 2001 ;[START_REF] Masini | Comparison of methodologies for determination of carboxylic and phenolic groups in humic acids[END_REF][START_REF] De Souza | Application of potentiometry to characterize acid and basic sites in humic substances: Testing the BEST7 program with a weak-acid mixture[END_REF].

TITRATIONS POTENTIOMETRIQUES

Les expériences de titrations ont été menées en cellule thermostatée à 25±0.2°C, en utilisant des solutions standardisées, sans CO 2 , de KOH (0.0999±0.0021M, à partir de KOH 0.5033M Sigma-Aldrich) et HNO 3 (0.1961±0.0002M, à partir de HNO 3 69% J.T. Baker), sous agitation et bullage d'azote saturé en eau et décarbonaté par barbotage dans KOH 1M. Les solutions standards de KOH et HNO 3 présentent une concentration en NaNO 3 de 0.10M. Le stand de micro-titrations utilisé est constitué de 2 Titrino 719 et d'un Titrino 716 Metrohm, reliés à une électrode de pH combinée (Mettler, Inlab422, référence: Ag/AgCl/KCl 3M), pilotés par le logiciel Tinet 2.4. L'électrode de pH est calibrée hebdomadairement par des solutions tampon HANNA 4.01, 7.01 et 10.01. 50mL de solutions de MO extraites sont utilisés pour les dosages acido-basiques. Ces solutions sont ajustées à pH 2 par ajouts de la solution standard d'HNO 3 , puis maintenues à ce pH pendant 20 minutes afin de les décarbonater. Ensuite, par des incréments variables de KOH standard (5 à 100µL), le pH est porté à 11.9 avec les conditions limites suivantes: 0.5mV.min -1 ou 120s. Une courbe de titration (de pH 2 à 12) est obtenue, contenant de 150 à 250 points (pH, V KOH ajouté). Seuls les points expérimentaux dont le pH est compris entre 2.5 et 11.5 sont utilisés pour la modélisation des propriétés d'acidité des MON, ceci afin de minimiser les éventuelles erreurs causées par la non-linéarité de l'électrode de pH pour les valeurs de pH extrêmes.

TITRATIONS VOLTAMMETRIQUES

Les concentrations initiales en métaux ont été déterminées, après digestion sous UV des solutions, par DPASV (Differential Pulse Anodic Stripping Voltammetry) par ajouts dosés. Ces valeurs sont nécessaires aux calculs permettant de calculer la spéciation. Pour permettre l'analyse d'une large gamme de concentrations en sites complexants présents sur les MON, les titrations des échantillons ont été réalisées en mode logarithmique (i.e ajouts de solutions standards à 1, 10, 100 et 1000µM en Cd ou Pb, par incréments constants en logarithme de la concentration en métal total). Pour cela, un montage de titration a été spécialement développé, basé sur un stand Metrohm 663VA, équipé de 4 burettes automatiques 765, relié à un PGSTATS12 (Ecochemie) et à un pH-mètre Metrohm 713 (électrode de pH Mettler Inlab422), contrôlé par un PC utilisant le logiciel GPES 4.9. La cellule de mesure, thermostatée à 25±0.2°C, est constituée d'une électrode de mercure à goutte pendante (SMDE), reliée à une électrode de référence Ag/AgCl/KCl 3M et à une électrode auxiliaire en graphite. Pratiquement, après 50 minutes sous agitation à 2000 rpm puis 10 minutes de purge sous N 2 saturé en eau et agitation, quatre gouttes de mercure sont éjectées, l'étape de réduction/pré-concentration (10 minutes de dépôts à -1.1V) se faisant sur la cinquième. Les titrations, effectuées en mode logarithmique, ont été dupliquées pour Cd et Pb afin d'estimer la reproductibilité. Chaque expérience dure approximativement 48 heures.

MODELISATION DES DONNEES EXPERIMENTALES

Le modèle de complexation utilisé dans cette étude, basé sur une complexation métal:ligand de type 1:1, repose sur l'utilisation d'une distribution discrète de ligands représentant les propriétés de complexation et d'acidité analysables des MON étudiées. Par opposition à la définition classique de la MON par des sites de complexation, chacun caractérisé par une seule concentration et une seule constante de complexation ou d'acidité, le concept de « quasi-particules » [START_REF] Sposito | Trace metals in contaminated waters[END_REF] a été utilisé. Une « quasi-particule » est une entité virtuelle, pas une molécule ou un site réel, définie par quelques propriétés comme des constantes de stabilité vis-à-vis de métaux traces et/ou du proton, des caractéristiques de fluorescence, .... facilement utilisable à une utilisation dans un modèle mathématique discret. Afin de modéliser les expériences réalisées, i.e. le comportement de la MON vis-à-vis des métaux traces et du proton, ainsi que les phénomènes de compétition métal/métal et métal/proton, trois types de « quasi-particules » ont été définis: I: sites spécifiques au proton (afin de tenir compte de la présence d'environ 100 fois plus de sites acides que de sites complexants vis-à-vis des métaux sur la MON [START_REF] Buffle | Complexation reactions in aquatic systems[END_REF]), définis par une constante d'acidité (K H I , C I ), II: sites échangeables vis-à-vis des métaux, prenant également en compte les compétitions métal/métal et métal/proton, définis par des constantes de stabilité vis-à-vis des métaux et du proton (K H II , K M,i II , C II pour i métaux) III: sites spécifiques, modélisant certaines interactions de la MON hautement spécifiques à un seul métal, sans compétition avec les autres métaux (K H III , K M III , C III ). L'ensemble de ces « quasi-particules » est appelé « chimio-type », celui-ci représente les propriétés des MON par un set de paramètres. [START_REF] Massart | Handbook of chemometrics and qualimetrics: Part A[END_REF], ajuste les valeurs des paramètres d'affinité sélectionnés ([L T ], K ML , K HL , …), donnant alors un nouveau set de paramètres au module de spéciation, ceci jusqu'à ce qu'une valeur d'erreur (bias) optimale soit atteinte (1000 à 10000 itérations). La figure 1 présente le principe de fonctionnement du logiciel PROSECE, subdivisé en deux parties explicitant respectivement les modules de calcul de spéciation et d'optimisation des paramètres de complexation et d'acidité des « quasi-particules » définies.

RESULTATS ET DISCUSSION

Propriétés d'acidité des MON

Les données des dosages acido-basique des AF et AH analysés ont été modélisées par un set de six « quasiparticules » de type I (ce nombre étant optimal). Trois de ces « quasi-particules » ont des constantes d'acidité proches de celles de sites carboxyliques (pK i = 4-5), les trois autres des valeurs proches de sites phénoliques (pK i = 9-10). Ce nombre et cette répartition de sites corroborent les résultats obtenus par [START_REF] Masini | Comparison of methodologies for determination of carboxylic and phenolic groups in humic acids[END_REF]. Avec une erreur moyenne inférieure à 0.5% pour des pH compris entre 3 et 11, le modèle mis en place a donc permis de reproduire correctement le comportement des MON étudiées vis-à-vis des protons. La figure 2

CONCLUSION

L'étude présentée a permis de mettre en évidence la possibilité de modéliser les propriétés d'acidité et de complexation de MON standards par un « chimio-type » ; propriétés analysées par des titrations acido-basiques et voltammétriques réalisées en mode logarithmique. Le « chimio-type », basé sur un modèle discret, est constitué d'un set de trois types de « quasi-particules ». Le premier type représente les sites acides présents sur la MON et n'est donc défini que par une concentration et une constante d'acidité. Le second décrit les sites échangeables, complexant les métaux traces, et prenant en compte les phénomènes de compétition métal/métal et métal/proton. Pour cela, les « quasi-particules » de type II sont donc définies par une concentration, des constantes de complexation vis-à-vis de tous les métaux traces étudiées et une constante d'acidité. Enfin un troisième type a été utilisé pour définir les sites spécifiques à un seul métal, sans compétition vis-à-vis des autres. L'utilisation du logiciel PROSECE permet l'optimisation de l'ensemble des propriétés de ces « quasi-particules » à partir des données expérimentales avec une erreur moyenne inférieur à 0.5% sur les pH et 1% sur les pM. Cette opération a été menée jusqu'au bout pour l'AFSR. Ce protocole expérimental et le modèle développé, appliqués à des échantillons naturels devraient permettre de modéliser le comportement de la MON et donc d'en apprécier le rôle dans le transport des métaux pour des systèmes complexes tels que des estuaires.

figure 1 :

 1 figure 1: représentation schématique du fonctionnement de PROSECE. figure 1: schematic representation of PROSECE functioning.

  figure 2 : Répartition des « quasi-particules » de type I (pK i et L iT ) optimisées pour les solutions d'AFSR, AHSR, AFLR et AHLR. figure 2: Optimised distribution of "quasi-particles" of type I (pK i and L iT ) for SRFA, SRHA, LRFA and LRHA solutions.

Figure 3 :

 3 Figure 3 : modélisation, par le « chimio-type » détaillé dans le tableau 1, des titrations logarithmiques de l'AFSR à pH 4.6 et 7.8 par Cd et Pb. Représentation log-log des concentrations mesurées (mes) et calculées (calc) en métal labile (pM) et erreurs (error) en fonction des concentrations en métal total (pM T ).Figure 3: Modelling of SRFA's logarithmic titrations, at pH 4.6 and 7.8 towards Cd and Pb, by a "chimio-type". Log-log plot of mesured (mes) and calculated (calc) labile concentrations of metals (pM) and error (error) function of total metal concentrations (pM T ).

Figure 3 :

 3 Figure 3 : modélisation, par le « chimio-type » détaillé dans le tableau 1, des titrations logarithmiques de l'AFSR à pH 4.6 et 7.8 par Cd et Pb. Représentation log-log des concentrations mesurées (mes) et calculées (calc) en métal labile (pM) et erreurs (error) en fonction des concentrations en métal total (pM T ).Figure 3: Modelling of SRFA's logarithmic titrations, at pH 4.6 and 7.8 towards Cd and Pb, by a "chimio-type". Log-log plot of mesured (mes) and calculated (calc) labile concentrations of metals (pM) and error (error) function of total metal concentrations (pM T ).

Tableau 1: "chimio-type", représentant le comportement de l'AFSR vis-à-vis des métaux traces et du proton, constitué de 6 "quasi-particules" de type I, 4 de type II et 2 de type III.


Cédric Garnier

Ivanka Pižeta
email: pizeta@irb.hr

Stéphane Mounier

Vlado Cuculi

Jean-Yves Yves Benaïm

Vlado Cuculić

An analysis of

Keywords: Voltammetry, Metal-ligand complexing parameters, Multiligand systems, Mono-and bi-ligand complexes, Mixed-ligand(s) complexes, PROSECE, DMS

ouverte pluridisciplinaire HAL, est destinée au dépôt

Introduction

As written in numerous papers, there has been a lot of work done and many efforts put in measuring metal-ligand interactions in model solutions and natural waters, while it is evident that the metals speciation is a good indicator of their bioavailability and hence, their toxicity [1][2][3][4][5][6][7]. Voltammetric methods with standard addition of metal ions of interest have shown to be one of the most non-destructive and subtle techniques, able to recognize metal-ligand complexation by directly measuring the free and labile fraction of metal and distinguishing free metal and labile from inert metal complexes [START_REF] Heyrovský | Principles of Polarography[END_REF][START_REF] Galus | Fundamentals of Electrochemical Analysis[END_REF][START_REF] Bond | Modern Polarographic Methods in Analytical Chemistry[END_REF][START_REF] Bard | Electrochemical Methods. Fundamentals and Applications[END_REF][START_REF] Town | Speciation Analysis by Electrochemical Methods, Handbook of Elemental Speciation: Techniques and Methodology[END_REF]. There are numerous models which deal with metal speciation with natural organic matter that have become sophisticated taking into account more and more components and parameters affecting the metal speciation and conditions when the bioavailable species could be present, as well as the models of natural organic matter description, which is heterogeneous and complex in its nature. They use large database knowledge about thermodynamic stability constants when discrete ligand approach is applied [START_REF] Lu | [END_REF][14][15][START_REF] Sposito | Chemical Equilibria and Kinetics in Soils[END_REF][START_REF] Tipping | [END_REF][18] or when continuous functions are associated with NOM properties and affinity to metal, major cations or proton [19][20][21][22][23], both needing a verification in experiments of various techniques. By the combination of such model-experiment approach, one can advance in characterization of samples of interest and obtain a more global speciation conception. In concrete applications, i.e. experimental verifications, the complexity of the system imposed the need for simplification [START_REF] Buffle | Analytical Chemistry[END_REF][START_REF] Dzombak | [END_REF], a stepwise approach. Indeed, it has been shown that many situations and experimental data could be explained and covered by a simplified conception of 1:1 metal-ligand complexes with their numerical indications of capacity and conditional complex stability constants [26][27][28][29][30][31][32][33]; on the other hand for a more complex conception there were no means of enough sensitive measurements and calculations, regarding voltammetric method approach [34,35]. When facing the complex structure of the real sample matrix, we always have to bear in mind that from single technique measurements (voltammetric or any other) one can get information about the so called conditional stability constants, which does not diminish the significance of the analyses, but necessitates further interpretation.

Ligands with strong complexing abilities will tend towards making bi-ligand complexes with trace metals, especially in the situations of low concentrations of the metals available [18,[START_REF] Buffle | Analytical Chemistry[END_REF]34,[36][37][38][39][40][41][42][43][44][45][46]. In natural water systems, which contain different ligands, mixed metal-ligand complexes are more probable to occur than the pure single metal-ligand ones [47]. Their quantitative characterization is not an easy task and is not often found in the literature. In contrast, theoretical calculations based on the knowledge of simpler systems, analyses of factors affecting the recognition of ternary complexes, and the reviews about particular ligand's behaviour in the environment (e.g. [48][49][50][51][52]) can be found. Mixed metal-ligand complexes were studied in order to improve the sensitivity of the voltammetry as electroanalytical method by multiple enhancement of metal-mixed ligands complexes adsorption onto mercury drop electrode and formation of metal-mixed ligands-surface complex, known as synergetic adsorption [53][54][55]. Also, this could be a possible type of process (among those well known) for metals removing from aqueous systems. Experimental results on synergetic adsorption phenomena of uranium and copper at the hanging mercury drop electrode have been published [53][54][55]. If we broaden the problematic to the formation of metal-ligand-surface complexes (surfaces found in natural waters and sediments and surface adsorptions acting as a ligand), there is a growing motivation for the study of ternary complexes [START_REF] Davis | Adsorption of Metal Complexes at Oxide and Related Surfaces[END_REF].

In our previous paper, a detailed analysis of the state-ofthe-art in metal complexing parameters determination in natural water systems for 1:1 metal-ligand systems has been given [START_REF] Garnier | [END_REF]. Starting from these results and using the optimal conditions described, possible to perform with our nowadays voltammetric experimental setup [58] (i.e. nonlinear, logarithmic standard additions) we wanted to analyze in more detail the ability of our experimental conditions and analyzing tools to distinguish the situations different from 1:1 metal-ligand situations that might occur in natural water systems, such as one and two 1:2 metal-ligand systems, as well as mixed ligand systems with the studied metal.

Simulation of experiment and formation of data sets for fitting

The simulation of experiment, i.e. numerical generation of data points that strictly imitate the distribution and the number of standard additions, that takes into account the existing sensitivity of the available instruments and noise distribution curve, is useful in predicting the recognition of information comprised in real, experimental data points, as shown previously [START_REF] Garnier | [END_REF]59]. According to the design of experimental standard additions, each point representing the sum of the free and the labile fraction of added metal was calculated by the program MINEQL [START_REF] Westall | MINEQL: a computer program for the calculation of chemical equilibrium composition of aqueous systems[END_REF].

The concentration ranges, the number of additions of titrated metal as well as the logarithmic type of additions were taken from the previous results [START_REF] Garnier | [END_REF], i.e. for one and two ligand cases which were analyzed in this work, the concentration range of added metal from 1 to 457 nmol L -1 is covered by 25 standard addition points with constant increments of pM T , i.e.

M T,i = 10 -(pM T,0 +(pM T,0 -pM T,n )×(i/n))
where n is 25 and M T,0 , M T,n and M T,i are the initial, the final and ith total metal concentrations, respectively. This design of titration, named logarithmic, has been the subject of the development of an automatic apparatus, with DPASV measurements of labile trace metal concentrations [START_REF] Garnier | [END_REF].

Considering the non-uniform repartition of the standard deviation of the measurements obtained on real logarithmic addition titration of cadmium [START_REF] Garnier | [END_REF], it seemed interesting to try reproducing that type of variation and apply it in this study. The noise added to the simulated data before fitting the best corresponds to:

e(M) = 1 exp((M/M ∞ ) × 10) -((e 0 -e ∞ -1)/(e 0 -e ∞ )) +e ∞
where M and M ∞ are the concentrations of metal added during and at the end of the titration, respectively, e 0 and e ∞ the noise limits at the beginning and end of the titration, 20 and 1%, respectively. So the added noise depends on metal addition and consists of adding random noise to the data points between -e(M) and +e(M). This process of data modification by noise addition is presented in Fig. 1.

Taking into account that analytical techniques used for the characterization of natural ligand complexing properties, especially when real samples are concerned, can seldom detect the concentration values lower than 0.1 nmol L -1 , this value was put as a limit to remove all the data points obtained by calculation which are lower than 0.1 nmol L -1 , i.e. the detection Fig. 1. An example of modification of data simulated by MINEQL after noise addition and removal of the points below 1 × 10 -10 mol L -1 . limit (DL) = 1 × 10 -10 mol L -1 . So, the number of points of one data set exposed to fitting was not always 25, it came down even to 8 for the experiments where the complexation intensity was the strongest. In Fig. 1, 14 points are left after noise application and detection limit restriction.

The existence of two ligands (L) was simulated, a strong (X) and a weak (Y) of 30 and 300 nmol L -1 concentrations, respectively. The position of the points selected for simulation in log(β 1 ) × log(β 2 ) space, representing the first and the second stability constants between metal (M) and ligand 2 ), is presented in Fig. 2. These values have been chosen taking into account the results obtained in the literature (or on real systems) [27][28][29][31][32][33][START_REF] Garnier | [END_REF][39][40][41][START_REF] Town | [END_REF]. The stability constant values chosen for simulation are summarised in Table 1.

(L = X or Y) (β 1 = [ML]/[M][L]; β 2 = [ML 2 ]/[M][L]
Regarding the formation constants of mixed complexes, i.e. MXY, it has been shown [50,51] that they should be related to the formation constants of the bi-ligand complexes by a relation of the following kind:

β XY = ω × (β 2X × β 2Y ) 1/2 ,
with ω ≥ 2. In our analysis, we used three different values of ω for the same combinations of formation constants as those used for simulation of single complexes, i.e. 2, 2 × 10 3 and 2 × 10 6 (Table 1). Four different groups of models were simulated using MINEQL: S1 = ML (for X and Y); S2 = ML + ML (for X and Y); S4 = MX + MX 2 + MY; and S5 = MX + MX + MY + MY 2 + MXY. The preset values of complexing parameters used for simulation can be read out from Fig. 2 and from Tables 2-5, respectively. For example: data set 7 stands for log β 1 = 10, log β 2 = 16, L = X = 30 nmol L -1 , the relevant results are found in Table 3.

For each group of simulated points, three forms of data sets were prepared, namely: (a) for all 25 originally calculated data points, (b) for data points after removal of the points lower than the selected determination limit DL, in our case DL = 1 × 10 -10 M (number of data points diminished up to 8), and (c) for data points obtained after removal of data <DL and application of the random noise.

Fitting strategy

Also resulting from our previous paper [START_REF] Garnier | [END_REF], the program PROSECE and its complete procedure was applied in all the fittings. Besides its proven suitability shown in the previous analyses, its structure is much more flexible and able to adapt itself to various fitting models used in this study compared to the other programs that could have been applied instead [59]. Its drawback is, however, the lack of fitted parameters confidence limits estimate.

In order to limit the total number of fittings, which in such analyses easily becomes too big and the obtained results not easy to elaborate, the schedule of fitting generally comprised the models similar or simpler than those used to simulate the data sets. Altogether five different fitting models were applied to data sets, namely, F1 = ML; F2 = ML + ML 2 ; F3 = MX + MY; F4 = MX + MX 2 + MY and F5 = MX + MX 2 + MY + MY 2 + MXY.

In the case of the "right" (matching) fitting model with the simulation model, one can discuss the errors of fitting, i.e. the differences between the preset and the obtained metal-ligand complexing parameters, but since in general we do not know which model is the right one, we should pay more attention to the other parameters of fitting, such as the amount and the shape of the residual of fitting. The parameter that also seemed to be interesting and drawing attention when analyzing the residuals was the number of zero crossings of the residual.

Results and discussion

Before we start the fitting procedure that should comprise several models, we should try to globally estimate the results by watching the shape of the data expressed in term of M f versus M T , pM f versus pM T and linearization modes of Scatchard [63] and/or Ružić-Van den Berg [35,[64][65][66], where M f stands for concentration of metal (free and labile) measured by voltammetry and M T for total metal added in standard addition procedure. That could give us the idea of the kind of fitting model to be used, or at least to discard some of the models. In the first presentation (M f versus M T ), one can try to estimate the total ligand capacity by noticing where the tangent to the data points intercepts X-axis, suiting the Chau method of linearization [67]. In the second presentation (pM f versus pM T ), if there is more than one inflection, and in the third presentation (linearization) if there is an indication of curvature of the data points, we can be sure that the structure is more complicated than 1:1 metal-ligand complexation.

The comparison of the three forms of data set results shows that difference is bigger between (a) and (b) than between (b) and (c). This is promising since the error function applied to the simulated data points is arbitrary in any case. Also, the detection limit of 10 -10 mol L -1 , used to remove the data points was rather severe; it could be expected to be lowered up to two orders of magnitude, depending on the analytical technique used. It has to be noticed that noise insertion was not that important for the final results concerning the goodness of fit.

As only the (c) case represents possible simulation of the experiment and as usually we do not know which model for fitting the experimental data is the right one, the most valuable information about the fit we have, besides logical parameter values, i.e. positive numbers for concentration and stability constants values, is the value of the residuals (res = sum(abs(ln(C calc )ln(C sim )))), which is the sum of absolute values of differences between the natural logarithms of calculated and simulated data points. The values of residuals are normalised to the number of data points in order to allow the comparison of the results for data sets of various lengths. It is interesting to compare consequent values obtained for the same set of simulated data by application of different fitting models.

The normalized residuals of all fittings for the data sets (c) (data points obtained after removal of data <DL and application of the random noise) of the four groups of simulated S1, S2, S4 and S5 models fitted by up to five F1-F5 fitting models, are shown in Fig. 3 and will be discussed separately. In general, at first sight, for the S1 and S2 cases the residuals are similar for all fittings, in which case, with the help of statistics we decide whether to accept a more complicated model or stay to the simpler one. For the S4 and S5 cases, we can notice that the application of the right model can in average be distinguished from the wrong ones by comparing the normalized residual values. In the tables the matching model is marked with an asterisk (*).

One-ligand systems (S1 and S2 simulations)

One-ligand system was presented with both ligands, L = X and Y, and their corresponding stability constants for monoligand as well as bi-ligand cases. Each of the simulated data sets was fitted by three different models, F1 = ML, F2 = ML + ML 2 and F3 = MX + MY.

S1 simulation

Six different data sets were considered, three for stronger (X) and three for weaker (Y) ligands; the preset values of complexing parameters, capacity and stability constants are given in Table 2 and can also be read out from Fig. 2. The results of (a) and (b) fittings with ML model were rather , d, e, f, g, h, i, j, k andl straightforward, almost perfectly matching the preset values, and are not surprising and to be discussed. The results of the simulated points fitting after removal of <DL points and noise insertion with three different models are shown in Table 2. Concerning normalized residuals (residuals divided by number of points), all the three fittings have similar residuals. In that case the simplest model should be accepted, which is the right decision according to the F-test [START_REF] Miller | Statistics and Chemometrics for Analytical Chemistry[END_REF]. It is interesting that for F2 model, the second constant β 2 has just been figured out by the fitting program, with almost no influence on the result of total ligand (except for the Sim. a case). The F3 model applied on strong ligands gives unacceptable results, but when applied on weaker ligands it just divides the total ligand in two groups of almost equal constants, which points that the decision of accepting the F1 model as the right one would be the correct choice.

S2 simulation

Eighteen different data sets were simulated, 9 for stronger (X) and 9 for weaker (Y) ligands. The results of data fitting from the (c) cluster with the same three models as in the S1 case are shown in Table 3. Looking only to the residuals provides no enough evidence to determine which model for fitting is the right one. On the contrary, for the set of stronger (X) ligands, lower residuals were obtained with a model offering two, rather than one ligand with two constants for mono-and bi-ligand complexes. In the group of S2 simulations, there is expressed the strongest impact of passing from fitting of all 25 data points (a), to fitting of data points without those eliminated by the DL (b) and fitting of noisy data (c). By analyzing only the (c) cluster, one can hardly decide about the right model and consequently about the values of the complexing parameters. In such a case the answer lies in redesigning the experiments either by measuring more replicates or refining and increasing the number of standard additions in order to increase the detection limit and lower the noise.

Two-ligand systems (S4 and S5 simulations)

Two-ligand systems were presented with the mixtures of ligands X and Y with the metal M and stability constants for mono-, bi-ligand cases as well as mixed metal-ligand complexes. Each of the simulated data sets was fitted by four or five various models, F1 = ML, F2 = ML + ML 2 , F3 = MX + MY, F4 = MX + MX 2 + MY and F5 = MX + MX 2 + MY + MY 2 + MXY.

S4 simulation

Fifteen different data sets were simulated following the S4 model, i.e. there is one stronger (X) ligand that forms mono-and bi-ligand-metal complexes and one weaker (Y) ligand that forms mono complex. According to their simulation names given in Table 4 and Fig. 2, one can read out their ligand concentrations and stability constants used. In Table 4, the obtained fitting parameters for the matching model as well as the values of the residuals for all the models applied are given together with the number of points of each data set for the (c) cluster.

In 9 out of the 15 cases the normalized residual was the smallest for the matching model. The weaker ligand concentration was retrieved with better score than the stronger one, with less than 10% of error, 11 out of the 15 and 3 out of the 15, respectively. Stability constant of the weaker ligand was well retrieved in all the 15 cases, while for the stronger ligand this happened in only 3 out of the 15 cases. Those positive cases were characterized at the same time with the number of points considerably high, i.e. the constants were the lowest of the proposed constant sets.

S5 simulation

Fifteen different data sets were simulated following the S5 model, i.e. there is one stronger ligand that forms monoand bi-ligand-metal complexes, one weaker ligand that forms mono-and bi-ligand-metal complexes, and there is one mixed-ligand complex. According to their simulation names given in Table 5, Fig. 2 and Table 1, one can read out their ligand concentrations and stability constants used. In Table the obtained fitting parameters for the matching model as well as the values of the residuals for all the models applied are given together with the number of points of each data set for the (c) cluster. In 10 out of the 15 cases the normalized residual was the smallest for the matching model, which is a better score than in the case of S4 simulation. However, the values of complexing parameters were estimated with higher tolerances. Ligand concentrations were underestimated, weaker ligands in higher extent than the stronger ones, in more than half of the cases. From five of the simulated stability constants the best retrieved constants after noising procedure were β 1Y and β 2Y , the constants for weaker mono-and bi-ligand-metal complexes. It is evident that for such rather complicated model, broader range of values of simulated parameters should be selected and the simulated data consequently treated in order to survey a specific rule and/or a range of more certain recognition of preset parameters.

Starting from these results, two characteristic cases were studied in more detail, i.e. the case d + 4, ω = 2 and the case f + 6, ω = 2 × 10 6 . Three-fold repetition of the experiments was simulated in order to test its influence on the correct retrieving of the parameters. The results given in Table brought to the conclusion that measurement repeating in great extent eliminates the influence of the noise (seen from the comparison of the corresponding (b) and (c) set of results), but cannot compensate for the missing of the data points <DL (comparison of a-type of data with b-and c-types in Table 6), which in fact carry the major part of information about the strong complex parameter values.

Data mining server rules

Handling with a big quantity of calculated data, we tried to possibly find a general rule that could help us recognize the right model of fitting, in our case simulated data, but in general, measured data of some unknown samples.

To that purpose we have prepared our data according to the rules of the data mining server (DMS) [START_REF] Gamberger | Data Mining Server[END_REF], where data analysis is performed based on knowledge induction by the inductive learning by logic minimization (ILLM) system. Using known data sets and their classification this program searches for classification rules that could be applied to the new data sets of the same kind. An attempt of mixing (a), (b) and (c) sets of data did not give any meaningful result, so we restricted ourselves to the (c) set. A table with 239 example rows (containing all the mentioned sets of simulated data fitted to all the mentioned models) was constructed, of them was found in the positive class, which in our case was a matching simulation and fitting model. Defined attributes were all preset parameter values of metal complexing capacity and their mutual differences and products, the number of data points in a particular set, the results of fitting including the amount of residual of fitting, and the residual normalized to the number of data points (res n ) as well as the number of zero crossings of the residual. A reasonable result is a rule obtained by the generalization parameter g = 100, which claims: fitting is good if res n < 2.53. This rule has a sensitivity of 100% and a specificity of 53.9%. Such a result was a confirmation of classical statistical approach to the fitting problem where the analysis of residuals is provided.

Another table (with 75 example rows, 15 out of which were positive, i.e. matching simulation and fitting model) was prepared with the examples of only S5, fitted with all the F1-F5 fitting models. A reasonable result is a rule obtained by the generalization parameter g = 20, which claims fitting is good if res n < 2.19 (sensitivity of 100% and specificity of 72.4%). When compared to the first table of simulated and fitting models all together, the second table with only one simulated model (the most complex in our analysis), shows more specificity for less generalization, which was expectable. Furthermore, the second table is a good model of a real experiment data treatment. By applying this rule, we could count on selecting all (100%) the matching models and 27.6% of the wrong models. Further combining with other disposable evidences such as critical logical inspection of the parameters values and error distribution (e.g. number of zero crossing of the residuals) could eventually bring us higher specificity.

Conclusion

The analyzed web of simulated experiments presuming and modelling the existence of more complex and probable bi-ligand and mixed metal-ligand complexes in aquatic systems, showed that with the nowadays state-of-the-art in measurement and data treatment techniques, it is possible to distinguish them from the simpler 1:1 metal-ligand systems in most of the analyzed cases.

This study showed that by means of voltammetric instrumentation supplied with automatic burettes and sophisticated data treatment that enables quick scanning through different fitting models, it is possible to conclude more precisely about the metal-ligand structure in a measured sample. However, for more accurate values of the complexing parameters, a more detailed study should be undertaken including repeating of experiments, and/or putting denser metal additions (in cases of strong ligands when larger number of first additions are so complexed to leave the free metal below the detection limit), considering that rather a subtle structure of dissolved metal-ligand complexes is tried to be distinguished.

Analyses by the method of data mining server confirmed the statistical methods of residuals analysis and have shown to be applicable in such simulation-fitting modelling systems as a prediction tool and help as one of the criteria in data interpretation. A more specific rule is to be expected with a denser web of modelled parameters, which could be the aim of some future work.

Fig. 2 .

 2 Fig. 2. The legend of distribution of metal first and second stability constants with stronger (X = 30 nmol L -1 ) and weaker (Y = 300 nmol L -1 ) ligand used to simulate the experiments.

Fig. 3 .

 3 Fig. 3. Residuals of fitting for data sets obtained after removal of data below the DL and application of random noise. S for simulation and F for fitting models: S1 = F1 = ML; S2 = F2 = ML + ML 2 ; F3 = MX + MY; S4 = F4 = MX + MX 2 + MY; S5 = F5 = MX + MX 2 + MY + MY 2 + MXY. Legends correspond to tables notation and are named from left to right for each F. S1: 1, 2, 3, a, b and c; S2: 4, 5, 6, 7, 8, 9, 10, 11, 12, d, e, f, g, h, i, j, k and l; S4: a + 4, a + 6, a + 8, a + 10, a + 12, b + 4, b + 6, b + 8, b + 10, b + 12, c + 4, c + 6, c + 8, c + 10 and c + 12; S5: d + 4 + ω = 2, d + 4 + ω = 2 × 10 3 , d + 4 + ω = 2 × 10 6 , f + 6 + ω = 2, f + 6 + ω = 2 × 10 3 , f + 6 + ω = 2 × 10 6 , h + 8 + ω = 2, h + 8 + ω = 2 × 10 3 , h + 8 + ω = 2 × 10 6 , j + 10 + ω = 2, j + 10 + ω = 2 × 10 3 , j + 10 + ω = 2 × 10 6 , l + 12 + ω = 2, l + 12 + ω = 2 × 10 3 and l + 12 + ω = 2 × 10 6 .

  Fig. 3. Residuals of fitting for data sets obtained after removal of data below the DL and application of random noise. S for simulation and F for fitting models: S1 = F1 = ML; S2 = F2 = ML + ML 2 ; F3 = MX + MY; S4 = F4 = MX + MX 2 + MY; S5 = F5 = MX + MX 2 + MY + MY 2 + MXY. Legends correspond to tables notation and are named from left to right for each F. S1: 1, 2, 3, a, b and c; S2: 4, 5, 6, 7, 8, 9, 10, 11, 12, d, e, f, g, h, i, j, k and l; S4: a + 4, a + 6, a + 8, a + 10, a + 12, b + 4, b + 6, b + 8, b + 10, b + 12, c + 4, c + 6, c + 8, c + 10 and c + 12; S5: d + 4 + ω = 2, d + 4 + ω = 2 × 10 3 , d + 4 + ω = 2 × 10 6 , f + 6 + ω = 2, f + 6 + ω = 2 × 10 3 , f + 6 + ω = 2 × 10 6 , h + 8 + ω = 2, h + 8 + ω = 2 × 10 3 , h + 8 + ω = 2 × 10 6 , j + 10 + ω = 2, j + 10 + ω = 2 × 10 3 , j + 10 + ω = 2 × 10 6 , l + 12 + ω = 2, l + 12 + ω = 2 × 10 3 and l + 12 + ω = 2 × 10 6 .

Table 1

 1 Values of stability constants of mono-, bi-and mixed complexes used for simulation in various combinations

	log(β 1X )	log(β 2X )	log(β 1Y )	log(β 2Y )	log(β XY )
	8	16	6	12	14.3, 17.3 and 20.3
	10	20	7	14	17.3, 20.3 and 23.3
	12	24	8	16	20.3, 23.3 and 26.3

Table 2

 2 Results obtained with PROSECE fitting of the six S1-type experiments using F1, F2 and F3 complexing models on c-type of data sets (n: number of simulated points after removal of points <DL, L T : ligand concentrations in nmol L -1 , res n : normalized residuals of PROSECE fitting)

	Ligand S1 sim. Preset values		F1 * = ML			F2 = ML + ML 2			F3 = MX + MY
			log(β 1 ) L T	n	log(β 1 ) L T	res n	log(β 1 ) log(β 2 ) L T	res n	log(β 1X ) X T	log(β 1Y ) Y T	res n
	X	1	8	3 0	2 5	8 .0	3 0 .0 0.951	8.0	1 2 .7	3 0 .2 0.952	8.0	30.0 4.2	245	0.949
		2	1 0	3 0	1 4	9 .9	3 0 .0 0.120 10.0	1 5 .0	3 0 .0 0.124 10.0	30.0 5.1	6.0 0.119
		3	1 2	3 0	1 1 1 2 .1	3 0 .0 0.173 11.5	1 7 .6	3 0 .0 0.173 12.7	30.0 -0.5	106	0.173
	Y	a	6	300	25	6.0	288	0.600	5.5	1 1 .7	550	0.598	6.2	19.7 6.0	271	0.600
		b	7	300	25	7.0	298	0.740	6.9	1 2 .4	319	0.740	7.2	19.7 7.0	280	0.740
		c	8	300	20	8.0	304	0.994	8.0	9 .4	304	0.994	8.0	18.5 8.0	285	0.994
	Table 3
	Results obtained with PROSECE fitting of the 18 S2-type experiments using F1, F2 and F3 complexing models on c-type of data sets; other parameters are the
	same as in Table 2
	S2 sim. Preset values			F1 = ML			F2 * = ML + ML 2			F3 = MX + MY
		log(β 1 ) log(β 2 ) L T	n	log(β 1 ) L T	res n	log(β 1 ) log(β 2 ) L T	res n	log(β 1X ) X T	log(β 1Y ) Y T	res n
	4	8	16	30	24	8.8	17	1.612	8.2	16.0	2 6 .6 0.848	8.8	6 .3	5.5	279	1.129
	5	8	20	30	14	10.8	15	0.185 10.8	17.8	1 4 .9 0.185 10.7	4 .9	2.4	292	0.185
	6	8	24	30	14	16.0	15	0.248 18.2	19.4	1 5 .0 0.248 14.7	5 .0 -63.0	1112	0.248
	7	10	16	30	14	10.0	30	0.453 10.0	14.5	3 0 .2 2.137 10.0	0 .0	3.0	0	0.453
	8	10	20	30	14	10.0	19	1.667 11.5	20.9	1 9 .7 1.483 10.3	7 .9	5.6	197	0.923
	9	10	24	30	14	16.5	15	0.326 19.5	11.3	1 5 .2 0.326 15.1	5 .2 -14.4	808	0.326
	10	12	16	30	11	10.7	30	0.157 10.6	9 .1	2 9 .8 0.156 11.1	9 .7	4.3	123	0.156
	11	12	20	30	14	16.4	30	0.123 20.1	8 .9	2 9 .9 0.123 11.2	9 .9	5.8	11	0.109
	12	12	24	30	14	10.0	20	1.631 11.2	20.7	2 0 .9 1.394 10.1	8 .7	5.6	200	0.908
	d	6	12	300	25	6.3	216	0.351	6.2	11.4	227	0.351	6.7	2 .1	4.2	202	0.353
	e	6	14	300	24	7.9	110	1.678	6.8	14.2	238	1.114	8.4	3 .4	5.9	203	1.171
	f	6	16	300	10	9.5	140	3.918	9.5	16.5	142	1.029	8.2	6 .5	5.7	206	1.470
	g	7	12	300	25	7.1	282	0.837	6.9	12.7	322	0.784	7.6	8 .9	5.1	199	0.777
	h	7	14	300	24	7.9	156	1.868	7.4	14.1	220	1.124	8.6	7 .9	6.1	188	1.513
	I	7	16	300	10	9.6	132	3.444	9.6	16.1	135	2.387	7.8	7 .3	5.3	217	0.944
	j	8	12	300	20	8.0	300	1.360	8.0	9 .6	300	1.360	8.0	4 .0	5.5	204	1.361
	k	8	14	300	19	8.1	276	1.623	2.8	13.8	755	1.277	8.8	1 .0	6.3	201	1.222
	l	8	16	300	9	9.4	158	5.398	9.2	16.4	184	2.373	6.6	0 .1	4.1	200	4.886

Table 4

 4 Results obtained with PROSECE fitting of the 15 S4-type experiments using the following fitting models: F1 = ML, F2 = ML + ML 2 , F3 = MX + MY and F4 = MX + MX 2 + MY; other parameters are the same as in Table2

	S4 sim. Preset values				F1, res n F2, res n F3, res n F4*
		log(β 1X ) log(β 2X ) X T log(β 1Y ) Y T	n			res n	log(β 1X ) log(β 2X ) X T	log(β 1Y ) Y T
	a + 4	8	16	30	6	300 24 3.716	2.233	0.914	0.732	8.5	.9	19.8 6.1	256
	a + 6	8	24	30	6	300 14 4.059	2.954	0.170	1.166 14.0	.5	21.6 5.7	475
	a + 8	10	20	30	6	300 14 5.877	5.833	0.419	0.428 10.2	.5	18.0 6.4	186
	a + 10	12	16	30	6	300 11 3.366	3.351	0.233	0.235 11.4	.1	30.4 6.0	302
	a + 12	12	24	30	6	300 14 5.797	5.758	0.410	0.419 10.3	.1	18.3 6.4	185
	b + 4	8	16	30	7	300 23 5.956	5.955	1.158	0.911	8.7	.9	16.8 7.0	317
	b + 6	8	24	30	7	300 14 6.910	6.910	0.337	0.844 13.7	.4	27.2 6.9	297
	b + 8	10	20	30	7	300 14 8.875	7.764	0.608	0.593 11.4	.3	15.9 7.0	302
	b + 10	12	16	30	7	300 11 5.767	5.767	0.917	0.516 11.2	.3	30.8 7.0	306
	b + 12	12	24	30	7	300 14 9.279	9.279	0.348	0.348 10.7	.8	17.6 7.0	317
	c + 4	8	16	30	8	300 19 1.231	1.231	1.062	1.126	8.0	.6	12.6 8.0	318
	c + 6	8	24	30	8	300 13 4.594	4.594	3.958	0.931 11.6	.4	16.0 8.0	304
	c + 8	10	20	30	8	300 13 4.592	4.592	4.313	0.743 15.1	.9	14.9 8.0	309
	c + 10	12	16	30	8	300 11 6.771	6.771	6.566	0.761 13.8	-2.1	30.4 8.0	299
	c + 12	12	24	30	8	300 13 4.993	4.993	4.693	1.061 10.4	.6	19.3 8.0

Table 5

 5 Results obtained with PROSECE fitting of the 15 S5-type experiments using the following fitting models: F1 = ML, F2 = ML + ML 2 , F3 = MX + MY, F4 = MX + MX 2 + MY and F5 = MX + MX 2 + MY + MY 2 + MXY; other parameters are the same as in Table2

	S5 sim. ω	n	F1, res n F2, res n F3, res n F4, res n F5 *
								res n	log(β 1X ) log(β 2X ) X T	log(β 1Y ) log(β 2Y ) log(β XY ) Y T
	d + 4	2	23	4.356	2.708	1.100	0.930	0.908 8.5	16.1	23.9 5.8	12.1	12.2
	d + 4	2 × 10 3 12	4.501	4.501	0.247	0.406	0.924 9.0	18.0	27.4 6.3	12.9	20.1
	d + 4	2 × 10 6 11	3.464	3.422	0.211	0.191	0.236 9.3	24.6	30.4 6.3	11.2	29.5
	f + 6	2	8	2.984	2.984	0.694	1.612	0.576 7.3	14.7	19.7 9.6	16.6	23.5
	f + 6	2 × 10 3	8	3.547	3.547	2.411	2.732	0.795 7.2	15.3	18.2 9.7	16.8	19.6
	f + 6	2 × 10 6	8	2.656	2.656	1.393	0.595	0.765 5.9	18.7	16.8 9.6	16.5	18.3
	h + 8	2	13	9.316	9.316	1.085	1.305	1.001 8.8	16.1	19.2 7.6	14.3	19.2
	h + 8	2 × 10 3 11	6.013	6.013	1.008	3.567	0.475 9.6	16.3	29.8 7.4	14.0	18.9
	h + 8	2 × 10 6 11	6.504	6.504	3.604	1.344	0.395 8.8	12.0	29.6 7.2	14.1	17.9
	j + 10	2	11	7.084	7.084	6.883	0.545	0.560 4.4	4.6	31.0 8.0	12.3	26.4
	j + 10	2 × 10 3 11	7.015	7.015	6.834	0.748	0.917 7.1	15.2	29.4 8.0	7.5	19.9
	j + 10	2 × 10 6 11	7.037	7.037	6.709	0.472	0.443 8.6	18.0	30.9 8.0	9.0	19.1
	l + 12	2	8	8.021	7.900	1.957	8.088	1.969 8.3	15.1	13.8 9.1	16.6	18.6
	l + 12	2 × 10 3	8	6.697	6.697	5.035	5.609	1.660 7.3	16.2	16.9 9.2	16.7	17.1
	l + 12	2 × 10 6	8	7.051	7.051	5.634	6.879	1.721 7.5	14.4	22.1 9.2	16.8	17.0

Table 6

 6 Results obtained with PROSECE fitting of two S5-type experiments using the fitting model F5; type: data type-see text; av: average of the three repetitions of the same simulated point with different random noise; other parameters are the same as in Table2

	Sim.	Fitting model	n	Type	n of meas.	res n	log(β 1X )	log(β 2X )	log(β 1Y )	log(β 2Y )	log(β XY )	X T	Y T
				Preset values			8.0	16.0	6.0	12.0	14.3	30.0	300
	d + 4ω = 2	F5 *		25	a	1	0.194	8.7	16.1	6.0	12.7	14.4	17.4	277
		F5 *		23	b	1	0.056	8.4	16.1	6.0	12.5	14.6	21.6	286
		F5 *		23	c	1	0.908	8.5	16.1	5.8	12.1	12.2	23.9	368
		F5 *		69	c	3	1.066	8.4	16.1	5.9	12.2	14.2	24.0	341
		F5 *		23 (av)	c	3	0.623	8.5	16.1	5.7	12.2	8.4	23.9	383
	Sim.		Fitting	n	Type	n of meas.	res n	log(β 1X )	log(β 2X )	log(β 1Y )	log(β 2Y )	log(β XY )	X T	Y T
			model	Preset values			8.0	24.0	6.0	16.0	26.3	30.0	300
	f + 6ω=2 × 10 6	F5 *	25	a	1	3.071 -3.6	25.2	5.4	16.1	26.3	31.2	264
			F5 *	8	b	1	0.583	6.7	16.6	9.6	16.7	24.2	17.5	156
			F5 *	8	c	1	0.765	5.9	18.7	9.6	16.5	18.3	16.8	152
			F5 *	24	c	3	1.326	7.1	15.9	9.6	16.7	20.9	19.9	157
			F5 *	8 (av)	c	3	0.663	7.2	16.1	9.6	16.7	21.6	19.2	157

Acknowledgement

The authors wish to thank the Ministry of Science, Education and Sport of the Republic of Croatia and the Ministry of National Education, Teaching and Research of France. Also we thank the Ph.D. school of the University of Sud Toulon Var for financial support of collaboration between France and Croatia and D. Gamberger for his help in explanation and application of DMS rules. We would like to thank two anonymous reviewers for their comments and suggestions.