Kinetic characterization of SPR-based biomarker assays enables quality control, calibration free measurements and robust optimization for clinical application

: Biomarker measurements are essential for the early diagnosis of complex diseases. However, many current biomarker assays lack sensitivity and multiplexing capacity, work in a narrow detection range and importantly lack real time quality control opportunities, which hampers clinical translation. In this paper, we demonstrate a toolbox to kinetically characterize a biomarker measurement assay using Surface Plasmon Resonance imaging (SPRi) with ample opportunities for real time quality control by exploiting quantitative descriptions of the various biomolecular interactions. We show an accurate prediction of SPRi measurements at both low and high concentrations of various analytes with deviations < 5% between actual measurements and predicted measurement. The biphasic binding sites model was accurate for fitting the experimental curves and enables optimal detection of heterophilic antibodies, cross-reactivity, spotting irregularities and/or other confounders. The toolbox can also be used to create a (simulated) calibration curve, enabling calibration-free measurements with good recovery, it allows for easy assay optimizations, and could help bridge the gap to bring new biomarker assays to the clinic.


Introduction
Early diagnosis of diseases can avoid progression, improve clinical outcomes and reduce healthcare costs [1].However, in the onset of diseases, changes between healthy subjects and patients are often minor and difficult to detect.Disease specific proteins, in this paper referred to as biomarkers, can give vital information on the disease progression from an early stage.In addition, the complexity of many diseases, and the heterogeneity of the patient populations, requires the measurement of several biomarkers simultaneous in order to acquire the desired differentiating capacity [2].Therefore, early disease diagnosis would benefit from sensitive, multiplex biomarkers assays.
The most important biomarker assay currently used in both research and clinic relies on the ELISA format or variations thereof.This tool has proven clinical effectivity and has shown to be reliable for over 30 years.However, it can measure only a single biomarker at a time and has a relatively small dynamic range, often requiring multiple dilution series [3].To be able to diagnose more complex diseases, several multiplex assays have been developed [2].These assays can be subdivided into 2D planar and bead suspension assays.Planar assays function similarly to ELISAs, with a sequential build up in complexity relying on one end-point measurement.The final read-out can be colorimetric, fluorescent, chemiluminescent or electro-chemiluminescent, with Searchlight (chemiluminescent) [4] and Mesoscale Discovery (electro-chemiluminescent) [5] as examples.Sensitivity of these assays is generally high, but the dynamic range is still limited [5,6].This is especially important for multiplex assays, as different biomarkers often exist in a wide variety of concentrations.Furthermore, these assays are shown to experience difficulty from spotting irregularities and unreliability, low precision [7,8], and lack of quality control opportunities [8,9].Suspension assays are other multiplex alternatives to the standard ELISA.These assays combine all reagents in suspension and therefore avoid the need for washing steps.Well-known examples are the flow cytometry bead arrays [6,10]and the Luminex™ assay [10].These assays show acceptable sensitivity and dynamic range [11].However, they suffer from large inter-assay variation at low concentrations [12] and inherent increase in cross-reactivity due to the suspension format [13].
Despite the development of several multiplex platforms, currently none are routinely used for patient diagnosis in a clinical setting [9].The main reason for this is the increasing complexity of the measurements as the number of biomarkers increases.For each biomarker added, additional cross-reactivity can occur, leading both to more difficult assay development and to more complex quality control [13].In addition, confounders in the sample matrix, for example heterophilic antibodies, can interfere in the interaction with any of the individual biomarkers [14].The black-box nature of the current multiplex assays, based on end-point measurements, only makes it extremely difficult to detect these confounders and hampers quality control, seriously limiting clinical use.
In order to advance the use of multiplex assays in clinical settings, more advanced methods for quality control with respect to end-point assays are needed.Potentially interesting in this respect are assays that enable real-time signal measurements.In combination with amplification steps, this can enable measurements in low fg/ml [15], with acceptable dynamic range [16] and multiplexing capabilities [17,18].We have recently described the development of a multiplex biomarker measurement system based on the real time surface plasmon resonance array imaging (SPRi) system.We have applied a sandwich antibody system followed with a neutravidin and biotinylated gold nanoparticle cascade to sequentially improve the detection signal [19].With this platform we have shown that we can quantitatively measure biomarkers in low fM (~fg/ml) concentration in combination with biomarkers up to 1 μg/ml in one assay depending on the quality of the capture antibody.In addition, this method proved that it is suited for multiplex measurements in a complex sample matrix such as serum and synovial fluid.This assay measures biomarker interactions in real-time allowing for extensive quality control and optimization opportunities for multiplex measurements.
In this paper, we describe the development of a simulation toolbox incorporating these opportunities.The method consists of a complete kinetic characterization of the biomarker cascade and the individual steps.First, we have fitted the individual cascade steps for test biomarkers IL-1β, IL-6, TNF-α and IFN-γ by applying both a simple 1:1 monophasic model kinetics model and a biphasic model of binding sites.We subsequently determined the relationship between the interaction steps over a concentration between 10 pg/ml and 10 ng/ml.This allowed us to fit the entire cascade and extrapolate the results accurately for all concentrations tested.We show that this toolbox enables improved quality control by simple detection of confounders, avoids the need of extensive calibration series and it can be used for rational assay optimization.Furthermore, besides saving time, it can reduce effort and costs.

Sensor preparation
Il-1β antibodies were immobilized on G-Type easy2spot sensors (Ssens, Enschede, The Netherlands) by reaction to free amines using the Wasatch microfluidics continuous flow spotter (Carterra, Salt Lake City, UT, US) for 30 min.We have used a gel-type sensor for the efficient use of the evanescent field and large binding capacity.The antibodies were prepared in a reaction buffer of 50 nM acetic acid (150 μl per spot).
Optimal reaction pH was 4.6 and assured concentration to negatively charged sensor and favoured coupling to primary amines.After successful spotting, the sensor was deactivated with 1% BSA in reaction buffer for 7 min and with 0.2 M ethanolamine at pH8.5 for another 7 min to reduce non-specific interactions.

Instrumentation
The SPRi measurements were performed on the IBIS MX96™ (IBIS Technologies, Enschede, the Netherlands).The instrument applies an angle-scanning method with automatic fitting to determine SPR shift with 0.5 Hz frequency.An automatic fluid-handling system controls a back-and-forth flow (20 μl, 30 μl/s) through a microfluidic flow cell (12 μl) to ensure minimal sample use.It can measure up to 96 spots simultaneously, making it highly suitable for multiplex experiments.The measurements were programmed using SUIT software (IBIS Technologies, Enschede, the Netherlands).The type of interaction, interaction times, samples, and regions of interest (ROIs) for the antibodies were set and a template was created that was loaded into IBIS data acquisition software.Before each experiment the angle offset was set to ensure wide dynamic detection range.SPRintX software was used for data collection and referencing.Data was subsequently exported to Matlab R2015A for evaluation using custom scripts which are available upon request.

SPRi enhancement cascade
SPRi signal enhancement cascades measurement was performed with IL-1β, IL-6, TNF-α and in a broad dynamic detection range as described before [19].Antibodies (cAb) were spotted at 5 μg/ml, 2.5 μg/ml, 1.25 μg/ml, 0.625 μg/ml, 0.3125 μg/ml and 0.15625 μg/ml, to achieve multiple ligand densities, with eight spots per concentration.Samples were dissolved in system buffer, containing PBS with 0.075% Tween80 and 0.5% BSA.Cytokines were measured at a concentration ranging from 100 fg/ml (~5 fM) to 1 μg/ml (~50 nM), spanning a dynamic detection range of 7 logarithms.The cascade interaction was performed as follows: First a cytokine (Cyt) was injected for 120 min, followed by a specific biotinylated detection antibody (dAb) at 5 μg/ml (33 nM) for 30 min, neutravidin (NeuAv) at 1.5 μg/ml (25 nM) for 15 min and a biotinylated gold nanoparticle (GNP) (40 nm diameter) at 77.69 mg/ml (0.2 nM) for another 15 min.After each interaction, the sensor was washed to reduce the non-specific signal and after each cascade, the sensor was regenerated using a double regeneration pulse for 30s, consisting of 200 mM phosphoric acid at pH 1.5.Every experiment starts with at least two blank measurements to obtain a steady baseline.In this cascade, each step sequentially grows the complex and therefore the SPRi signal.This is summarized in table S1.

Fitting of cascade steps
The cascade steps were kinetically characterized using a biphasic model (a 1:1 monophasic model, figure s1, and mass transport limitation model, figure s2, did not explain the curvature), accounting for two binding sites, one with relatively high and one with relatively low affinity.For this we used exponential equation (1) [20,21] representing the biphasic binding sites model: In this equation R max1 and R max2 are the binding capacities (RU), k a1 and k a2 the association rate constants (M − 1 s − 1 ), and k d1 and k d2 the dissociation rate constants (s − 1 ) of the ligand species 1 and 2, respectively, c is the analyte concentration (M) and t is time from start of the interaction (s).
We have used custom Matlab scripts to fit the different models (available upon request).For the 1:1 and biphasic model, we have used the "lsqcurvefit" function embedded in Matlab.As input, we have provided time, R signal, and concentrations.Additionally, we applied upper and lower limits to the fits, determined by initial guess and physiological feasibility.We applied different fitting strategies for the analyte interaction step and the cascade enhancement steps.We fitted the antibodyanalyte kinetics globally over the analyte concentrations of 1 ng/ml, 10 ng/ml and 100 ng/ml and over the spotting concentrations range (5 μg/ ml-0.15μg/ml) simultaneously.k a (1,2) and k d (1,2) were shared variables, c a concentration dependent input, and R max (1,2) was allowed to vary between the different ligand spotting concentrations.
The enhancement cascade steps interact with the previously captured ligandanalyte complex.This complex can be considered as a ligand with a specific binding capacity (ligands defined in table s1).The variable in these interactions is now not the 'analyte' concentration, but the ligand density.For these steps we applied global fitting over all signals from the analyte concentration range from 10 pg/ml to 10 ng/ml and the spotting concentrations range (5 μg/ml -0.15 μg/ml) simultaneously (30 signals in total per step).k a (1,2) and k d (1,2) were shared variables, c a cascade step dependent input, and R max (1,2) was allowed to vary between the different ligand spotting concentrations.
The simple 1:1 fits and fits accounting for biphasic behavior resulted in a large number of R max (1,2) values for each cascade step depending on the ligand density.It is likely that binding is dominated by high affinity spots at low ligand densities leading to a high Rmax1/Rmax2 ratio [21].This would resemble a simple 1:1 fit.At higher binding capacities the number of lower affinity spots could increase, partially due to steric effects, and could lead to a low Rmax1/Rmax2 ratio.We analyzed these relationships in detail in supplemental figures s6-s8.

Relationship between cascade steps
The relationship between the cascade steps was determined by first acquiring R sat values (RU end -RU start ) at an analyte concentration range from 1 pg/ml to 10 ng/ml over a spotting concentration range from 5 μg/ml -0.15 μg/ml.R sat values for the analyte step at low concentrations (below 1 ng/ml) were below measurable threshold.For this specific subset, simulated R sat values from the biphasic model were used.We subsequently determined the relationship of the R sat values between the analyte-detection antibody, detection antibody -neutravidin, and neutravidin -gold nanoparticle steps.

Prediction of enhancement cascade
We applied the information on the kinetics of the individual cascade steps, the correlation between R max and R sat and the relationship between the cascade steps to simulate the enhancement cascade.To achieve this, we first simulate the analyte step for three low [1 pg/ml, 3.3 pg/ml, 10 pg/ml] and three high [100 pg/ml, 1 ng/ml and 10 ng/ml] IL-1β concentrations at a spotting concentration of 5 μg/ml by either the 1:1 model or the biphasic model (specific analyte concentrations for IL-6, TNF-α and IFN-γ are described in supplemental figure s12-s14).The resulting R sat,analyte is used to calculate R sat,dAb by applying the fitted model.R max , dAb(1,2) is derived from the R sat , dAb by the appropriate model and the detection antibody step can be simulated.R sat,dAb is subsequently used to calculate R sat,NeuAv , which is used to derive R max,NeuAv (1,2) .The same steps are again used to calculate the parameters for the GNP step.This leads to the complete simulation of the enhancement cascade for the concentrations tested for the four biomarkers of interest.

Calibration free measurements
The accurate prediction of the enhancement cascade allows for calibration free measurements.In principle only knowledge of the R max of the analyte binding to the ligand spot is required to predict the entire cascade.To demonstrate the potential, we have fitted the curves at a specific step in the enhancement cascade (Cyt, dAb, NeuAv or GNP) and reverse calculated the concentration.We have performed this for a concentration range from 1 pg/ml to 10 ng/ml at spotting concentration of 5 μg/ml by either the 1:1 model or the biphasic model.The procedure is as follows: For example, a fit is performed at the GNP step resulting in R max (1,2) , this can be recalculated to R sat,GNP , which can be recalculated to R sat,NeuAv , to R sat,dAb and to R sat,analyte .From this R sat,analyte the concentration of the analyte can be calculated based on the kinetics data.Using this method, we have created a simulated calibration curve.We have compared the simulated calibration curve to the concentrations calculated from calibration curve derived via the standard method in biomarker assays (fit according to equation ( 1)) and to the actual concentrations.

Assay optimization
The kinetic characterization of the enhancement cascade allows for optimization experiments using simulations.We have demonstrated this by adjusting the concentrations of the individual cascade steps and their association times.This results in enhancement cascades with different timings and signals.We have subsequently optimized the cascade to obtain maximal signal (by adjusting concentration) and optimized for both time and signal (by adjusting concentration and association times).

Kinetic characterization of biomarker assay
We have used the SPRi enhancement cascade that was described before [19], as a model biomarker assay that can be fully kinetically characterized.For our model reactions, we have measured IL-1β, IL-6, TNF-α and IFN-γ in a concentration range between 100 fg/ml and 1 μg/ml over six spot densities.A representative measurement for IL-1β at 10 ng/ml at a spotting concentration of 5 μg/ml is shown in Fig. 1A.Fig. 1A shows the sequential increase in SPRi signal in the cascade, from 25 RU after analyte interaction to 125 RU after detection antibody, 475 RU after neutravidin and to 4500 RU after gold nanoparticle interaction.This results in a total signal increase of 5, 15 and 180 times respectively.This large signal increase enables a lower limit of detection of 50 fg/ml in a dynamic range of more than 7 logarithms or > 10,000,000 times.Each step of the enhancement cascade is governed by specific kinetics and can as such be kinetically characterized.We applied a biphasic binding sites model (Fig. 1B and S2, S5) to characterize the antibodyanalyte interaction step (1:1 Langmuir interaction model shown in S1, S4 and S5).The figure shows that the antibody-antigen interaction can be fitted with small residuals with the two binding sites model (residuals increase slightly at high concentration, 100 ng/ ml, yet stay below 5%).This may indicate the presence of two different affinity sites, but more likely it is expected to indicate a range of affinities caused by surface heterogeneities.The fits provided a K D1 of 3.7 nM and a K D2 of 7.6 nM.The weighted K D was 6.9 nM which closely resembles the K D of the 1:1 fit (K D = 5.1 nM) and the 1:1 fit K D determined with extensive single cycle kinetics (SCK) experiment (K D = 1.7 nM).In the Supplemental Fig. S4, the fits for the other spot densities are shown.Similarly, the enhancement steps (1C-E) could be fitted with small residuals (more detailed results in figure s5).We have tested a more complex model with three ligand species, but this did not result in further reduced residuals (data not shown).Therefore, we suggested that these residuals might have a different origin.
We looked further into the two distinct binding species in each cascade step: one with faster association and one with slower association kinetics.We have assessed the correlation between the R max (1, 2) of the ligand species and total ligand density over a set of IL-1β and spotting concentrations.Here, we found a clear non-linear relationship between the R sat signal (and thus ligand density) and the R max of the ligand species.In all cascade steps we found a slow transition from a dominance of one ligand species at low total density to the second ligand species at higher total density (more details in supplemental figures s6-s8).
After the characterization of the individual cascade steps, the relationship between these steps was determined (Fig. 1F-H).The 4-parameter logistic model best described the relationship between the cascade steps, both the curvature at higher concentrations (missed by linear model, Fig. S9) and lower concentrations (missed by the exponential model, Fig. S10).The curvature is especially prominent in the relationship between analyte vs detection antibody and neutravidin vs gold nanoparticle and is characterized by a declining slope, indicating saturation.

Prediction of enhancement cascade
The characterization of the individual enhancement steps and their inter relationships allows us to predict the enhancement cascade at a specific IL-1β concentration and spotting concentration(Fig.2, prediction 1:1 model in Fig. S11).
Fig. 2 shows that the enhancement cascade can be predicted with good accuracy for both very low (deviation 2.9%) and high (deviation 4.3%) concentrations by applying the parameters from the biphasic binding sites model and the relationship between the cascade steps.This is in contrast with the 1:1 model that could only predict the correct curvature for the low concentrations (Fig. S11).These accurate predictions allow for extensive quality control opportunities.For example, spotting irregularities can be detected when deviation from expected curvature occurs at the capture antibody-analyte interaction.Similarly, deviation in the detection antibody steps suggest possible crossreactivity, where absence of curvature might indicate possible heterophilic antibodies.These confounders are easily missed in end-point assays and these quality control possibilities are a strength of our approach.
To show the universal application of this method, we additionally performed cascade characterizations of IL-6, TNF-α and IFN-γ.We demonstrate that these cascades can be accurately predicted using the same fitting procedures, and cascade step correlations (see figure s12-s14).In figure s15 we show how deviations in fitting parameters influence the accuracy of the cascade predictions.To show reproducibility, we repeated the affinity measurement 3 times in identical setting.We show that the fitting is reproducible with CV<15% (Table S6).

Calibration free measurements
In addition to extensive quality control, the accurate prediction of the enhancement cascade also allows for calibration free measurements.This can be achieved by reversing the steps required to predict the cascade at a specific IL-1β concentration.A fit at any step in the cascade can then be used to acquire the IL-1β concentration.The results are shown in Fig. 3 (results for 1:1 fit in Fig. S16).
Fig. 3 shows that the IL-1β concentration can be calculated from the simulations accurately until a lower limit of detection (LLoD) that is dependent on the step in the enhancement cascade (LLoD are 600 pg/ml, 200 pg/ml, 75 pg/ml and 1 pg/ml for analyte, detection antibody, neutravidin and gold nanoparticle step respectively), below this LLoD the experimental calibration curve is more accurate.However, the simulated calibration curve is as good as the normal calibration curve to determine concentrations in the linear range (mean recovery, defined as measured concentration divided by standard concentration x 100%, in linear range for gold nanoparticle is 108% in simulated calibration curve vs 94% for standard curve).This avoids the need for extensive calibration measurements, saving time and money.

Assay optimization
The kinetic characterization of the enhancement cascade allows us to easily optimize the biomarker assay without requiring large number of experiments.We can simply adjust concentrations and times of individual steps and simulate the cascade.As a proof of concept, we have performed two assay optimizations, one optimizing the signal and one optimizing both time and signal.The results are shown in Fig. 4.
In Fig. 4A the simulations to increase the total signal in the enhancement cascade is shown.Here, we explored the parameter space by increasing the concentration of each component of the enhancement cascade by 5 times (detection antibody, neutravidin and gold nanoparticle).This leads to an improvement in signal by ~25% (signal increase from 125 to 160 RU) and thus a higher sensitivity.Similarly, we optimized the enhancement cascade for both time and signal (Fig. 4B).When a 5 times higher concentration was used, this resulted in a reduction of assay time by 50% (180 min-90 min) with no signal decrease, improving the applicability for point of care applications.To further optimize the assay for the intended application, we simulated the limit of detection based on assay time (Fig. 4c).Based on this curve optimal assay times can easily be chosen, short (suitable for point of care, but with reduced LLOD), or long (suitable for standard diagnosis, with improved LLOD).This proof of concept demonstrates that optimizations for biomarker assays can be easily designed using our method, which can then be subsequently experimentally validated.

Discussion
Current biomarker assays lack required sensitivity, dynamic range, multiplexing capacity [2] and essential quality control making them unsuitable for large scale clinical applications [7].We have previously Fig. 3. Predictions of enhancement cascade can be used for calibration free measurements.A-D) IL-1β calculated concentration is plotted against the actual IL-1β concentration in the biomarker, detection antibody, neutravidin and gold nanoparticle step respectively.The dashed line shows the actual concentration, x shows calculated concentration based on the experimental calibration curve and the squares the calculated concentration based on the simulations (Results for 1:1 langmuir model are shown in Fig. S16).When signals fall below the calibration curves, calculated concentrations are set to 0. Signals above the calibration curve can also not be measured and are not shown.Additionally, due to fitting procedure calibration curves are only accurate between input data points (10 pg/ml -100000 pg/ml).Together, this means the calibration curves are limited between the LLoD and the HLoD, which are dependent on cascade step.. developed a biomarker assay based on SPRi [12] that enables biomarker measurements in multiplex with both high sensitivity and specificity in a broad dynamic detection range.
While kinetic modelling of individual SPR interactions has been extensively studied [22], in this paper we demonstrate that combining 1) full characterization of individual cascade steps and 2) the relationship between these steps in biomarker assays, enables extensive quality control, calibration-free measurements and allows for simple optimizations.The same steps can be applied to any biomarker assays using real time technologies, both label free and with enhancement tags.This can help bridge the gap and finally bring these technologies to the clinic.
Kinetic characterization of the biomarker assay showed a biphasic behaviour was required to accurately fit the interactions at both high and low concentrations.The 1:1 monophasic model strongly deviated at high spotting concentration and high analyte concentration, and interestingly no mass transport limitation was observed, despite high k a (>10 5 M − 1 s-1 ) in dAb, NeuAv and GNP steps [14].
These results indicate that there are distinct binding species in all steps of the enhancement cascade (both with high and with low affinity [13]).At low concentrations only the high affinity species will be occupied resulting in 1:1 behaviour, while at higher concentrations both high and low affinity species will become occupied leading to a biphasic curvature.This behavior was confirmed by the data showing that the R max of a binding species was highly correlated to the ligand (or cascade complex) density.The heterogeneity in the enhancement cascade can be caused by several processes.The covalent immobilization process of the capture antibody is random and can lead to differences in orientation and steric freedom leading to different affinities [23].This works through the cascade and partly explains the heterogeneity of subsequent steps.However, other factors unique to the architecture of each individual step could play a role in the enhancement steps.The biotinylation of the detection antibody is achieved by similar random chemistry, potentially leading to different species with specific steric biotin availability and antibody affinity [24].The interaction between the neutravidin and the gold nanoparticle is also limited by steric effects, especially at high concentrations, further increasing heterogeneity.This diversity would suggest a continuum of binding affinities instead of just two, as was described by Ref. [14].However, tests using three or four distinct species did not result in better fits for any of the steps.Therefore, because of the good fits and by applying Occam's razor, we decided to continue to use the biphasic binding model with two binding species.
Similar biphasic behaviour and heterogeneity can be expected in other biomarker assays to which our method can be applied.ELISAs for example use physical immobilization of antibodies, in which random orientation can be expected.Additionally, random biotinylation of detection antibodies will also play a role.For biomarker assays using SPR platforms, chemiluminescence, electro-luminescence, electrochemical or bead suspension methods similar or additional heterogeneity can be expected.Therefore, our approach is likely necessary and suitable for characterization of a broader range of assay types.
The relationship between the enhancement cascade steps could be captured accurately with the 4-parameter logistic model (and moderately with an exponential model), indicating a linear relationship between the cascade steps at relatively low concentrations that reaches a saturation at higher concentrations.This might be related to the saturation of the higher affinity binding spots, altering binding kinetics.Despite this, the good fits demonstrate that the cascade was well defined over the entire concentration range that was tested.As similar fits are standard use in biomarker assays, it is likely that these cascade relationships are universally applicable [25].
The kinetic characterization allowed accurate simulations of the enhancement cascade at both low and high IL-1β concentrations with average deviations of less than 5%.Additionally, accurate predictions were shown for IL-6, TNF-α (with small deviations due rebinding effect) and IFN-γ demonstrating general functionality of a broad range of antibodies and affinities.Accurate fitting of curvature using biphasic model (in contrast to 1:1 Langmuir fitting) enables detection of confounders that affect binding characteristics.This enables extensive quality control opportunities.Cross-reactivity is an important limitation in all multiplex assays, caused by cross-reactivity of capture-and detection antibodies to off-target biomarkers or other components of the sample matrix.This will result in a lower or higher signal, which will be interpreted as a lower or higher biomarker concentration.This is notoriously difficult to detect [15] and can lead to false positives or negatives that can be devastating in a clinical setting [10].Using our method, we can detect this cross-reactivity through knowledge of the kinetics.When a capture-or detection antibody non-specifically cross-reacts with a biomarker, this will impact the different binding curves of both the capture antibodybiomarker and the biomarkerdetection antibody interactions [26,27].In addition, the relationship between the enhancement cascade steps will start to deviate from expectation.Similarly, common confounders such as heterophilic antibodies and spotting irregularities will translate into experimentally determined binding curves which deviate from predicted curves.These deviations can be automatically detected leading to accurate quality control both for assay development and clinical application.In future work these applications will be experimentally verified.
The cascade simulations were additionally used to create a simulated calibration curve, resulting in accurate concentration measurements (108% recovery) vs a standard calibration curve (94% recovery).This indicates that the kinetic characterization of the biomarker assays allows for calibration free measurements.The only requirement is knowledge of the binding capacity of the capture antibody spot.This knowledge can be gained by a simple triple biomarker injection requiring very limited time.According to Mehand et al. [28] it might even be possible to achieve calibration free measurements with this method, without sensor specific knowledge on the ligand density, by determining this density in real time.This reduced need for calibrators can save time and reagents and can be very useful, especially in point of care applications.
Finally, we showed in proof of concept that kinetic characterization enables simple offline assay optimization.Resulting in signal improvement of 25% or assay time reduction of 50%.These simulations need to be tested in real experiments, but show that simple changes in concentration or assay time can be easily simulated enabling rational assay optimization rather than trial and error.This can save a lot of time, reagents and money in optimizing biomarker assays, especially when multiple markers are added and assays become more complex.
Currently, multiplex biomarker measurements fail to reach routine clinical application.This is caused in part by increased complexity when measuring multiple, and especially a large number of biomarkers.At high multiplexing, confounders, such as cross-reactivity, spotting irregularities and heterophilic antibodies are increasingly difficult to detect in end-point assays.
We have demonstrated that the kinetic characterization of optical biomarker assays can lead to extensive quality control opportunities, can lead to calibration free measurements, and can help in assay optimization.We have shown this for the enhancement cascade that we have developed, but the principles are applicable to all real-time assay technologies.The method works well for complex kinetics with a wide range of affinities (nM for antibodies and fM for neutravidin/GNP) and can be applied in label free and other nanoparticle tags assays.In addition, when knowledge on the specific kinetics of antibodybiomarker pairs is gathered using a real time technology, this can also easily be applied to optimize traditional end-point assays, such as the standard ELISA and multiplex assay.
We will make our scripts publicly available, so researchers can apply our method to their assay of choice.In principle, the scripts can be applied to automatically detect confounders for routine clinical use.This can help optimize these assays and enable useful quality control possibilities.In turn this can hopefully lead to widespread use of these assays in a clinical setting where they are desperately needed.

Conclusion
Biomarker assays are essential for the early diagnosis of diseases, but are currently lacking for clinical application.In this paper, we have developed a toolbox to kinetically define biomarker assays and show that this has large potential for quality control, calibration free measurements and assay optimizations.To achieve this characterization, we have performed three steps.1) Kinetically defining the individual assay steps of the amplification cascade using a biphasic binding sites model, 2) Defining the correlation of ligand binding species to interaction signal and 3) determining the relationship between cascade steps.Combined this resulted in an accurate prediction of SPRi measurements at both low and high concentrations of analytes with deviations <5% between actual measurements and predicted measurement.Deviation between prediction and real time measurements could point to cross-reactivity, heterophilic antibodies, spotting irregularities and/or other confounders.This method is applicable to all real-time assay technologies and can also be applied to optimize standard end-points assays.Finally, we show in a proof of concept that these simulations can be used to easily optimize signal (25% increase) while reducing assay time (>50% reduction), compared to the original assay, without requiring extensive experimental work.The theoretical framework of our assay can be applied to any real time biomarker assay relying on molecular interactions.It can also be used to optimize standard end-point assays.Therefore, this could be a useful tool to bring the desired biomarker assays for early diagnosis one step closer to the clinic.

Funding sources
Dutch Arthritis Association and Dutch Research Council.

Fig. 1 .
Fig. 1.A: The sequential build up in SPRi signal is shown in the enhancement cascade for IL-1β.Det.Ab = Detection antibody, NeuAv = Neutravidin, GNP = Gold nanoparticle.B) The fits for antibody-antigen interaction are shown over the concentrations 1, 10 and 100 ng/ml at a spotting concentration of 5 μg/ml (~1900 RU).The residuals, i.e. the difference between fit and actual measurements are shown in the subplot below.C-E) the fits of the enhancement cascade steps are shown.Fitting parameters are shown in table s3.F-H) Relation between the enhancement cascade steps.From left to right the relationship between analyte vs detection antibody, detection antibody vs neutravidin and neutravidin vs gold nanoparticle is shown.The signals for specific IL-1β concentration over the six spot densities is grouped by colour.The fits depict the results from the 4-parameter logistic model.(For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)

Fig. 2 .
Fig. 2. Prediction of the enhancement cascade.The predictions (dashed lines) are shown in comparison with the SPRi signals (solid lines) for three high (A-B) and three low (C-D) IL-1β concentrations (spotting concentration: 5 μg/ml).The analytedetection antibody step (A and C) and neutravidingold nanoparticle step (B and D) are shown separately to increase clarity.(For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)

Fig. 4 .
Fig. 4. Simulations to increase assay signal or reduce assay time for the enhancement cascade.A) The simulation to improve signal is shown compared to the standard cascade at an IL-1β concentration of 100 pg/ml.B) At the same IL-1β concentration the optimization for both time and signal is shown.C) The relationship between limit of detection and assay time is shown.